Data

Morphological data

Morphological measurements (eye diameter and body length) were taken from a mixture of fresh and preserved specimens. In total, 454 sergestid specimens were sampled from 15 species belonging to 11 genera and two subgroups (Sergestes and Sergia).

# Full dataset for ontogenetic allometry ------

#import data
specimens_raw <- data.frame(read.csv("../Data/specimen_data.csv", header=TRUE, na.strings=c("", "NA", " ", stringsAsFactors = FALSE)))

#tidy data
specimens <- specimens_raw %>%
  mutate(Genus = gsub('\\s+', '', Genus)) %>%
  mutate(Family = str_to_title(Family)) %>%
  mutate(Subgroup = str_to_title(Subgroup)) %>%
  mutate(Genus = str_to_title(Genus)) %>%
  mutate(genus_species = as.factor(paste(Genus, Species, sep = "_"))) %>%
  mutate(pres_bin = recode_factor(Preservation, "ethanol" = "preserved", "paraformaldehyde" = "preserved", "fresh" = "fresh"))


# Species means for evolutionary allometry ------

#import data
species_raw <- data.frame(read.csv("../Data/species_data.csv", header=TRUE, na.strings=c("", "NA", " ", stringsAsFactors = FALSE)))

#tidy data
species <- species_raw %>%
  mutate(Genus = gsub('\\s+', '', Genus)) %>%
  mutate(genus_species = as.factor(paste(Genus, Species, sep = "_"))) 

# Show sampling -----

#compile number of shrimp sampled by species
counts <-ddply(specimens, .(specimens$Subgroup, specimens$Genus, specimens$Species), nrow) 

#create table column names
names(counts) <- c("Subgroup", "Genus","Species","Sampled") 

#generate scrolling table in RMarkdown
kable(counts[ , c("Subgroup", "Genus","Species","Sampled")], caption = "Species and sampling effort for morphological data from family Sergestidae.") %>% 
                      kable_styling(full_width = F) %>% 
                      collapse_rows(columns = 1, valign = "top") %>%
                      scroll_box(height = "400px") 
Species and sampling effort for morphological data from family Sergestidae.
Subgroup Genus Species Sampled
Sergestes Allosergestes pectinatus 14
Allosergestes sargassi 38
Deosergestes corniculum 18
Deosergestes henseni 58
Neosergestes edwardsii 7
Parasergestes armatus 30
Parasergestes vigilax 17
Sergestes atlanticus 15
Sergia Challengerosergia hansjacobi 43
Challengerosergia talismani 32
Gardinerosergia splendens 42
Phorcosergia grandis 51
Robustasergia robusta 23
Robustosergia regalis 38
Sergia tenuiremis 28

Molecular phylogeny

Below is the sergestid molecular phylogeny generated by Charles Golightly that we pruned to the 15 species in the dataset for comparative analyses.

#Import phylogeny
tree_orig <- read.tree(file = "../Data/sertestid_tree_rooted")

#Make copy of tree to manipulate
tree_edit <- tree_orig

#Remove locality tags from phylo tip labels 
n <- 3 #which _ to end pattern at
pat <- paste0('^([^_]+(?:_[^_]+){',n-1,'}).*') #pattern
tree_edit$tip.label <- sub(pat, '\\1', tree_edit$tip.label) #remove tags from tips

#Make list of taxa to drop (in tree but not in dataset)
drops <- setdiff(tree_edit$tip.label, as.character(species$Binomial))

#Drop unwanted tips from phylogeny
tree.pruned <- drop.tip(phy = tree_edit, tip = drops) 

#Check and reorder tree
is.rooted(tree.pruned) 
is.binary.tree(tree.pruned) 
tree <- ladderize(tree.pruned) 

#Make row names of the species the phylogeny tip labels
rownames(species) <- species$Binomial

#Reorder species data to match phylogeny order
species <-ReorderData(tree, species, taxa.names="row names")

#check that phylogeny tips and data match exactly (if they match will return "OK")
name.check(tree, species)

#make tree ultrametric
tree <- chronopl(tree, 1)
#plot tree of sampled species with original tip labels
plot.phylo(tree, #phylogeny to plot
           type = "phylogram", #how to visualize phylogeny
           show.tip.label = TRUE, #whether to show tip labels/species names
           cex = 1.0, #text size
           no.margin = TRUE, 
           use.edge.length = TRUE,
           edge.width = 2, #thickness of phylogeny branches
           label.offset = 0.1) #how far away from tips to put labels

Ontogenetic eye-body allometry in Sergestidae

Eye-body allometry in each species

We first examine intraspecific/ontogenetic eye diameter-body length allometry across 15 species of deep-sea sergestid shrimp. Allometric relationships were fit using standardised major axis regression via the smatr package. For each species, standard model checks, model fits, and an interactive plot of data and fits are shown. Hover over the plots to see individual data values (this is super helpful for checking weird outliers that may be typos in the source datasheets).

We can also use SMA regressions to test whether allometry estimated from fresh vs. preserved specimens yield significantly different scaling relationships. However, since the sample sizes within species are not always very large, sometimes this results in very low sample sizes for each treatment, which can be problematic. For now, I have juts colored by preservation type and we can discuss how we want to handle these going forward.

col_pres <- c("ethanol" = "#a6cee3",
              "fresh" = "#b2df8a",
              "paraformaldehyde" = "#1f78b4")

Allosergestes pectinatus

# Subset data -----
apec <- specimens %>% filter(genus_species == "Allosergestes_pectinatus")

# Fit SMA model ------
sma_apec <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = apec, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_apec, which = "default",type = "o") 
plot(sma_apec, which = "residual",type = "o")
plot(sma_apec, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_apec)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = apec, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##             elevation    slope
## estimate    -2.908342 2.009753
## lower limit -3.816721 1.427844
## upper limit -1.999963 2.828814
## 
## H0 : variables uncorrelated
## R-squared : 0.6929183 
## P-value : 0.00022059
#save coefficients of fits as object
cc_sma_apec <- data.frame(coef(sma_apec))

# Make plot ------
plot_apec <- ggplot(apec, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Allosergestes pectinatus") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_apec, aes(intercept = cc_sma_apec[1,1], slope = cc_sma_apec[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_apec[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_apec[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_apec)

Allosergestes sargassi

Note: this species has a big outlier that we are thinking we should remove. Measured onboard ship, no way to check data

# Subset data -----
asar <- specimens %>% filter(genus_species == "Allosergestes_sargassi")

# Fit SMA model ------
sma_asar <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = asar, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_asar, which = "default",type = "o") 
plot(sma_asar, which = "residual",type = "o")
plot(sma_asar, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_asar)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = asar, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##             elevation    slope
## estimate    -2.212583 1.474440
## lower limit -2.806943 1.104000
## upper limit -1.618222 1.969178
## 
## H0 : variables uncorrelated
## R-squared : 0.246605 
## P-value : 0.0015179
#save coefficients of fits as object
cc_sma_asar <- data.frame(coef(sma_asar))

# Make plot ------
plot_asar <- ggplot(asar, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Allosergestes sargassi") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_asar, aes(intercept = cc_sma_asar[1,1], slope = cc_sma_asar[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_asar[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_asar[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_asar)

Challengerosergia hansjacobi

# Subset data -----
chan <- specimens %>% filter(genus_species == "Challengerosergia_hansjacobi")

# Fit SMA model ------
sma_chan <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = chan, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_chan, which = "default",type = "o") 
plot(sma_chan, which = "residual",type = "o")
plot(sma_chan, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_chan)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = chan, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##             elevation     slope
## estimate    -1.516521 1.0399481
## lower limit -1.879267 0.8284462
## upper limit -1.153776 1.3054464
## 
## H0 : variables uncorrelated
## R-squared : 0.4712721 
## P-value : 3.731e-07
#save coefficients of fits as object
cc_sma_chan <- data.frame(coef(sma_chan))

# Make plot ------
plot_chan <- ggplot(chan, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Challengerosergia hansjacobi") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_chan, aes(intercept = cc_sma_chan[1,1], slope = cc_sma_chan[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_chan[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_chan[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_chan)

Challengerosergia talismani

# Subset data -----
ctal <- specimens %>% filter(genus_species == "Challengerosergia_talismani")

# Fit SMA model ------
sma_ctal <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = ctal, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_ctal, which = "default",type = "o") 
plot(sma_ctal, which = "residual",type = "o")
plot(sma_ctal, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_ctal)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = ctal, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##             elevation     slope
## estimate    -1.564958 1.0316438
## lower limit -1.815606 0.8835032
## upper limit -1.314310 1.2046238
## 
## H0 : variables uncorrelated
## R-squared : 0.825775 
## P-value : 6.5338e-13
#save coefficients of fits as object
cc_sma_ctal <- data.frame(coef(sma_ctal))

# Make plot ------
plot_ctal <- ggplot(ctal, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Challengerosergia talismani") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_ctal, aes(intercept = cc_sma_ctal[1,1], slope = cc_sma_ctal[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_ctal[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_ctal[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_ctal)

Deosergestes corniculum

# Subset data -----
dcor <- specimens %>% filter(genus_species == "Deosergestes_corniculum")

# Fit SMA model ------
sma_dcor <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = dcor, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_dcor, which = "default",type = "o") 
plot(sma_dcor, which = "residual",type = "o")
plot(sma_dcor, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_dcor)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = dcor, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##             elevation     slope
## estimate    -1.337193 0.8105777
## lower limit -1.618827 0.6598417
## upper limit -1.055559 0.9957482
## 
## H0 : variables uncorrelated
## R-squared : 0.8471468 
## P-value : 6.2959e-08
#save coefficients of fits as object
cc_sma_dcor <- data.frame(coef(sma_dcor))

# Make plot ------
plot_dcor <- ggplot(dcor, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Deosergestes corniculum") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_dcor, aes(intercept = cc_sma_dcor[1,1], slope = cc_sma_dcor[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_dcor[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_dcor[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_dcor)

Deosergestes henseni

# Subset data -----
dhen <- specimens %>% filter(genus_species == "Deosergestes_henseni")

# Fit SMA model ------
sma_dhen <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = dhen, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_dhen, which = "default",type = "o") 
plot(sma_dhen, which = "residual",type = "o")
plot(sma_dhen, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_dhen)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = dhen, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##             elevation     slope
## estimate    -1.327868 0.8376931
## lower limit -1.616552 0.6763970
## upper limit -1.039183 1.0374523
## 
## H0 : variables uncorrelated
## R-squared : 0.3519041 
## P-value : 9.2412e-07
#save coefficients of fits as object
cc_sma_dhen <- data.frame(coef(sma_dhen))

# Make plot ------
plot_dhen <- ggplot(dhen, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Deosergestes henseni") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_dhen, aes(intercept = cc_sma_dhen[1,1], slope = cc_sma_dhen[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_dhen[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_dhen[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_dhen)

Gardinerosergia splendens

Note: big outlier

# Subset data -----
gspl <- specimens %>% filter(genus_species == "Gardinerosergia_splendens")

# Fit SMA model ------
sma_gspl <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = gspl, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_gspl, which = "default",type = "o") 
plot(sma_gspl, which = "residual",type = "o")
plot(sma_gspl, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_gspl)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = gspl, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##              elevation     slope
## estimate    -1.2559422 0.8778027
## lower limit -1.6211316 0.6687602
## upper limit -0.8907528 1.1521882
## 
## H0 : variables uncorrelated
## R-squared : 0.2574859 
## P-value : 0.00060377
#save coefficients of fits as object
cc_sma_gspl <- data.frame(coef(sma_gspl))

# Make plot ------
plot_gspl <- ggplot(gspl, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Gardinerosergia splendens") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_gspl, aes(intercept = cc_sma_gspl[1,1], slope = cc_sma_gspl[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_gspl[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_gspl[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_gspl)

Neosergestes edwardsii

Note: there is a big outlier in this species where the eye size is actually bigger than the body length (seems impossible) - flagging to check. Also, sample size is low and narrow range of body sizes so this species will probably not make the cut for developmental allometry but is great for the evolutionary comparisons

# Subset data -----
nedw <- specimens %>% filter(genus_species == "Neosergestes_edwardsii")

# Fit SMA model ------
sma_nedw <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = nedw, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_nedw, which = "default",type = "o") 
plot(sma_nedw, which = "residual",type = "o")
plot(sma_nedw, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_nedw)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = nedw, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##               elevation      slope
## estimate    -0.02725658 -0.1654676
## lower limit -0.19419020 -0.3524820
## upper limit  0.13967704 -0.0776764
## 
## H0 : variables uncorrelated
## R-squared : 0.4782384 
## P-value : 0.085234
#save coefficients of fits as object
cc_sma_nedw <- data.frame(coef(sma_nedw))

# Make plot ------
plot_nedw <- ggplot(nedw, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Neosergestes edwardsii") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_nedw, aes(intercept = cc_sma_nedw[1,1], slope = cc_sma_nedw[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_nedw[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_nedw[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_nedw)

Parasergestes armatus

# Subset data -----
parm <- specimens %>% filter(genus_species == "Parasergestes_armatus")

# Fit SMA model ------
sma_parm <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = parm, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_parm, which = "default",type = "o") 
plot(sma_parm, which = "residual",type = "o")
plot(sma_parm, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_parm)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = parm, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##             elevation     slope
## estimate    -1.660526 1.0318360
## lower limit -2.051439 0.8006561
## upper limit -1.269613 1.3297664
## 
## H0 : variables uncorrelated
## R-squared : 0.5613314 
## P-value : 1.9015e-06
#save coefficients of fits as object
cc_sma_parm <- data.frame(coef(sma_parm))

# Make plot ------
plot_parm <- ggplot(parm, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Parasergestes armatus") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_parm, aes(intercept = cc_sma_parm[1,1], slope = cc_sma_parm[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_parm[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_parm[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_parm)

Parasergestes vigilax

# Subset data -----
pvig <- specimens %>% filter(genus_species == "Parasergestes_vigilax")

# Fit SMA model ------
sma_pvig <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = pvig, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_pvig, which = "default",type = "o") 
plot(sma_pvig, which = "residual",type = "o")
plot(sma_pvig, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_pvig)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = pvig, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##             elevation    slope
## estimate    -3.256619 2.328849
## lower limit -4.339295 1.631234
## upper limit -2.173943 3.324807
## 
## H0 : variables uncorrelated
## R-squared : 0.5634777 
## P-value : 0.00051661
#save coefficients of fits as object
cc_sma_pvig <- data.frame(coef(sma_pvig))

# Make plot ------
plot_pvig <- ggplot(pvig, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Parasergestes vigilax") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_pvig, aes(intercept = cc_sma_pvig[1,1], slope = cc_sma_pvig[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_pvig[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_pvig[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_pvig)

Phorcosergia grandis

# Subset data -----
pgra <- specimens %>% filter(genus_species == "Phorcosergia_grandis")

# Fit SMA model ------
sma_pgra <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = pgra, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_pgra, which = "default",type = "o") 
plot(sma_pgra, which = "residual",type = "o")
plot(sma_pgra, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_pgra)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = pgra, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##             elevation     slope
## estimate    -1.463079 0.9333209
## lower limit -1.625863 0.8440017
## upper limit -1.300296 1.0320925
## 
## H0 : variables uncorrelated
## R-squared : 0.8768029 
## P-value : < 2.22e-16
#save coefficients of fits as object
cc_sma_pgra <- data.frame(coef(sma_pgra))

# Make plot ------
plot_pgra <- ggplot(pgra, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Phorcosergia grandis") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_pgra, aes(intercept = cc_sma_pgra[1,1], slope = cc_sma_pgra[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_pgra[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_pgra[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_pgra)

Robustasergia robusta

# Subset data -----
rrob <- specimens %>% filter(genus_species == "Robustasergia_robusta")

# Fit SMA model ------
sma_rrob <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = rrob, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_rrob, which = "default",type = "o") 
plot(sma_rrob, which = "residual",type = "o")
plot(sma_rrob, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_rrob)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = rrob, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##              elevation     slope
## estimate    -1.2367827 0.8547738
## lower limit -1.6029678 0.6709735
## upper limit -0.8705977 1.0889227
## 
## H0 : variables uncorrelated
## R-squared : 0.7097726 
## P-value : 4.5859e-07
#save coefficients of fits as object
cc_sma_rrob <- data.frame(coef(sma_rrob))

# Make plot ------
plot_rrob <- ggplot(rrob, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Robustasergia robusta") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_rrob, aes(intercept = cc_sma_rrob[1,1], slope = cc_sma_rrob[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_rrob[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_rrob[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_rrob)

Robustosergia regalis

# Subset data -----
rreg <- specimens %>% filter(genus_species == "Robustosergia_regalis")

# Fit SMA model ------
sma_rreg <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = rreg, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_rreg, which = "default",type = "o") 
plot(sma_rreg, which = "residual",type = "o")
plot(sma_rreg, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_rreg)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = rreg, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##             elevation     slope
## estimate    -1.317150 0.8937637
## lower limit -1.526402 0.7755937
## upper limit -1.107898 1.0299381
## 
## H0 : variables uncorrelated
## R-squared : 0.8227993 
## P-value : 4.2961e-15
#save coefficients of fits as object
cc_sma_rreg <- data.frame(coef(sma_rreg))

# Make plot ------
plot_rreg <- ggplot(rreg, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Robustosergia regalis") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_rreg, aes(intercept = cc_sma_rreg[1,1], slope = cc_sma_rreg[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_rreg[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_rreg[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_rreg)

Sergestes atlanticus

# Subset data -----
satl <- specimens %>% filter(genus_species == "Sergestes_atlanticus")

# Fit SMA model ------
sma_satl <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = satl, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_satl, which = "default",type = "o") 
plot(sma_satl, which = "residual",type = "o")
plot(sma_satl, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_satl)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = satl, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##             elevation    slope
## estimate    -2.629861 1.805583
## lower limit -3.474390 1.301400
## upper limit -1.785333 2.505095
## 
## H0 : variables uncorrelated
## R-squared : 0.6905247 
## P-value : 0.00012413
#save coefficients of fits as object
cc_sma_satl <- data.frame(coef(sma_satl))

# Make plot ------
plot_satl <- ggplot(satl, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Sergestes atlanticus") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_satl, aes(intercept = cc_sma_satl[1,1], slope = cc_sma_satl[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_satl[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_satl[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_satl)

Sergia tenuiremis

# Subset data -----
sten <- specimens %>% filter(genus_species == "Sergia_tenuiremis")

# Fit SMA model ------
sma_sten <- sma(formula = Eye_Diameter ~ Body_Length, 
            data = sten, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            alpha = 0.05)

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_sten, which = "default",type = "o") 
plot(sma_sten, which = "residual",type = "o")
plot(sma_sten, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_sten)
## Call: sma(formula = Eye_Diameter ~ Body_Length, data = sten, log = "xy", 
##     method = "SMA", alpha = 0.05) 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Coefficients:
##              elevation     slope
## estimate    -1.1996991 0.7854599
## lower limit -1.6293146 0.5673032
## upper limit -0.7700836 1.0875090
## 
## H0 : variables uncorrelated
## R-squared : 0.3252106 
## P-value : 0.0015327
#save coefficients of fits as object
cc_sma_sten <- data.frame(coef(sma_sten))

# Make plot ------
plot_sten <- ggplot(sten, aes(x = Body_Length, y = Eye_Diameter, col = Preservation)) + 
  geom_point(size = 2, alpha = 0.8) + 
  scale_color_manual(values = col_pres, name = "Sample state", breaks = c("ethanol", "fresh", "paraformaldehyde")) +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  ggtitle("Sergia tenuiremis") + 
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), plot.title = element_text(face = "italic")) +
  geom_abline(data = cc_sma_sten, aes(intercept = cc_sma_sten[1,1], slope = cc_sma_sten[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_sma_sten[2,1], digits = 2), sep = " "), size = 3, col = "black") +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_sma_sten[1,1], digits = 2), sep = " "), size = 3, col = "black")

ggplotly(plot_sten)

Figure for developmental eye-body scaling across sergestids

#make figure panels
fig.a <- plot_parm + theme(legend.position = "none", axis.title.x = element_blank())
fig.b <- plot_pvig + theme(legend.position = "none", axis.title.x = element_blank(), axis.title.y = element_blank())
fig.c <- plot_dcor + theme(legend.position = "none", axis.title.x = element_blank(), axis.title.y = element_blank())
fig.d <- plot_dhen + theme(legend.position = "none", axis.title.x = element_blank())
fig.e <- plot_asar + theme(legend.position = "none", axis.title.x = element_blank(), axis.title.y = element_blank())
fig.f <- plot_apec + theme(legend.position = "none", axis.title.x = element_blank(), axis.title.y = element_blank())
fig.g <- plot_chan + theme(legend.position = "none", axis.title.x = element_blank())
fig.h <- plot_ctal + theme(legend.position = "none", axis.title.x = element_blank(), axis.title.y = element_blank())
fig.i <- plot_pgra + theme(legend.position = "none", axis.title.x = element_blank(), axis.title.y = element_blank())
fig.j <- plot_sten + theme(legend.position = "none", axis.title.x = element_blank())
fig.k <- plot_rrob + theme(legend.position = "none", axis.title.x = element_blank(), axis.title.y = element_blank())
fig.l <- plot_rreg + theme(legend.position = "none", axis.title.x = element_blank(), axis.title.y = element_blank())
fig.m <- plot_gspl +theme(legend.position = "none")
fig.n <- plot_satl + theme(legend.position = "none", axis.title.y = element_blank())
fig.o <- plot_nedw + theme(legend.position = "none", axis.title.y = element_blank())

# arrange panels in figure
plots <- plot_grid(fig.a, fig.b, fig.c, fig.d, fig.e, fig.f, fig.g, fig.h, fig.i, fig.j, fig.k, fig.l, fig.m, fig.n, fig.o, #list of plots to arrange in grid
           align = 'vh', #align horizontally and vertically
           axis = 'lb',
           labels = c("A", "B", "C", "D", "E", "F", "G", "H", "I","J", "K", "L", "M", "N", "O"), #panel labels for figure
           hjust = -1, #adjustment for panel labels
           nrow = 5) #number of rows in grids

#extract legend
leg <- get_legend(fig.a + theme(legend.position="bottom"))


#view figure
plot_grid(plots, leg, nrow = 2, rel_heights = c(1, .01), align = 'vh')

#export figure
png("../Figures/dev-allometry.png", width = 10, height = 12, units = "in", res = 500)
plot_grid(plots, leg, nrow = 2, rel_heights = c(1, .05), align = 'vh')
dev.off()

The plots here are scaled to fit the data, so they all look similar but the slopes are actually quite variable and I think pretty interesting. This could be worth exploring for a chunk of the paper, and we might consider examining variation in developmental scaling slopes as a factor of phylogeny and things like whether each species exhibits ontogenetic descent.

Test whether species differ in ontogenetic allometry

The smatr package also allows you to test whether groups exhibit significant differences in allometric scaling (slopes or intercepts). Below, we again use standardized major axis regression to model the scaling of eye diameter with body length across all specimens measured and test whether we see significant differences in developmental scaling across species. This is a bit redundant with above, but will provide statistical evidence of the (very big) difference in slopes among species.

*Note: As we only have 7 specimens of Neosergestes edwardsii and one is an odd outlier while the other 6 are all very similar in body length, we don’t have the power to fit developmental eye-body scaling for that species. Thus, I’ve removed it from the specimen dataset before fitting this model.

This test shows that sergestid species significantly differ in their allometric scaling relationships between eye diameter and body length (p < 0.0001). It also yields a table of pairwise comparisons so you can see exactly which species differ in slope and which do not.

#Filter/prep dataframe
specimens_sma <- specimens %>%
  filter(genus_species != "Neosergestes_edwardsii") %>% # Remove N. edwardsii from specimen data
  mutate(species_text = as.factor(paste(Genus, Species, sep = " "))) # Add labels for species 

# Model for developmental scaling across specimens grouped by species
sma_species <- sma(formula = Eye_Diameter ~ Body_Length * genus_species, 
            data = specimens_sma, 
            log = "xy", #sets both x and y variables as logged
            method="SMA", #defines SMA as model
            multcomp = TRUE, # adds pairwise comparisons between groups
            multcompmethod = "adjusted", # adjusts p-value for multiple comparisons
            alpha = 0.05) 

#plot fit, residuals, qq
par(mfrow = c(2,2)) #make plot window 2x2
plot(sma_species, which = "default",type = "o") 
plot(sma_species,which = "residual",type = "o")
plot(sma_species, which = "qq", type = "o")
par(mfrow = c(1,1)) 

#view fits
summary(sma_species)
## Call: sma(formula = Eye_Diameter ~ Body_Length * genus_species, data = specimens_sma, 
##     log = "xy", method = "SMA", alpha = 0.05, multcomp = TRUE, 
##     multcompmethod = "adjusted") 
## 
## Fit using Standardized Major Axis 
## 
## These variables were log-transformed before fitting: xy 
## 
## Confidence intervals (CI) are at 95%
## 
## ------------------------------------------------------------
## Results of comparing lines among groups.
## 
## H0 : slopes are equal.
## Likelihood ratio statistic : 60.07 with 13 degrees of freedom
## P-value : 5.0941e-08 
## ------------------------------------------------------------
## 
## Results of multiple comparisons among groups.
## 
## Test for pair-wise difference in slope :
##                 genus_species_1              genus_species_2      Pval
## 1      Allosergestes_pectinatus       Allosergestes_sargassi 0.9999999
## 2      Allosergestes_pectinatus Challengerosergia_hansjacobi 0.1953505
## 3      Allosergestes_pectinatus  Challengerosergia_talismani 0.1066269
## 4      Allosergestes_pectinatus      Deosergestes_corniculum 0.0066684
## 5      Allosergestes_pectinatus         Deosergestes_henseni 0.0112631
## 6      Allosergestes_pectinatus    Gardinerosergia_splendens 0.0371468
## 7      Allosergestes_pectinatus        Parasergestes_armatus 0.2210213
## 8      Allosergestes_pectinatus        Parasergestes_vigilax 1.0000000
## 9      Allosergestes_pectinatus         Phorcosergia_grandis 0.0223771
## 10     Allosergestes_pectinatus        Robustasergia_robusta 0.0190940
## 11     Allosergestes_pectinatus        Robustosergia_regalis 0.0154883
## 12     Allosergestes_pectinatus         Sergestes_atlanticus 1.0000000
## 13     Allosergestes_pectinatus            Sergia_tenuiremis 0.0174772
## 14       Allosergestes_sargassi Challengerosergia_hansjacobi 0.9967671
## 15       Allosergestes_sargassi  Challengerosergia_talismani 0.9484813
## 16       Allosergestes_sargassi      Deosergestes_corniculum 0.0892585
## 17       Allosergestes_sargassi         Deosergestes_henseni 0.1796345
## 18       Allosergestes_sargassi    Gardinerosergia_splendens 0.6120191
## 19       Allosergestes_sargassi        Parasergestes_armatus 0.9980406
## 20       Allosergestes_sargassi        Parasergestes_vigilax 0.9881695
## 21       Allosergestes_sargassi         Phorcosergia_grandis 0.2909931
## 22       Allosergestes_sargassi        Robustasergia_robusta 0.3429079
## 23       Allosergestes_sargassi        Robustosergia_regalis 0.2046918
## 24       Allosergestes_sargassi         Sergestes_atlanticus 1.0000000
## 25       Allosergestes_sargassi            Sergia_tenuiremis 0.3456937
## 26 Challengerosergia_hansjacobi  Challengerosergia_talismani 1.0000000
## 27 Challengerosergia_hansjacobi      Deosergestes_corniculum 0.9999497
## 28 Challengerosergia_hansjacobi         Deosergestes_henseni 1.0000000
## 29 Challengerosergia_hansjacobi    Gardinerosergia_splendens 1.0000000
## 30 Challengerosergia_hansjacobi        Parasergestes_armatus 1.0000000
## 31 Challengerosergia_hansjacobi        Parasergestes_vigilax 0.0332808
## 32 Challengerosergia_hansjacobi         Phorcosergia_grandis 1.0000000
## 33 Challengerosergia_hansjacobi        Robustasergia_robusta 1.0000000
## 34 Challengerosergia_hansjacobi        Robustosergia_regalis 1.0000000
## 35 Challengerosergia_hansjacobi         Sergestes_atlanticus 0.4856371
## 36 Challengerosergia_hansjacobi            Sergia_tenuiremis 0.9999999
## 37  Challengerosergia_talismani      Deosergestes_corniculum 0.9966572
## 38  Challengerosergia_talismani         Deosergestes_henseni 0.9999891
## 39  Challengerosergia_talismani    Gardinerosergia_splendens 1.0000000
## 40  Challengerosergia_talismani        Parasergestes_armatus 1.0000000
## 41  Challengerosergia_talismani        Parasergestes_vigilax 0.0148260
## 42  Challengerosergia_talismani         Phorcosergia_grandis 1.0000000
## 43  Challengerosergia_talismani        Robustasergia_robusta 1.0000000
## 44  Challengerosergia_talismani        Robustosergia_regalis 1.0000000
## 45  Challengerosergia_talismani         Sergestes_atlanticus 0.2845038
## 46  Challengerosergia_talismani            Sergia_tenuiremis 0.9999977
## 47      Deosergestes_corniculum         Deosergestes_henseni 1.0000000
## 48      Deosergestes_corniculum    Gardinerosergia_splendens 1.0000000
## 49      Deosergestes_corniculum        Parasergestes_armatus 0.9999985
## 50      Deosergestes_corniculum        Parasergestes_vigilax 0.0006718
## 51      Deosergestes_corniculum         Phorcosergia_grandis 1.0000000
## 52      Deosergestes_corniculum        Robustasergia_robusta 1.0000000
## 53      Deosergestes_corniculum        Robustosergia_regalis 1.0000000
## 54      Deosergestes_corniculum         Sergestes_atlanticus 0.0168318
## 55      Deosergestes_corniculum            Sergia_tenuiremis 1.0000000
## 56         Deosergestes_henseni    Gardinerosergia_splendens 1.0000000
## 57         Deosergestes_henseni        Parasergestes_armatus 1.0000000
## 58         Deosergestes_henseni        Parasergestes_vigilax 0.0012337
## 59         Deosergestes_henseni         Phorcosergia_grandis 1.0000000
## 60         Deosergestes_henseni        Robustasergia_robusta 1.0000000
## 61         Deosergestes_henseni        Robustosergia_regalis 1.0000000
## 62         Deosergestes_henseni         Sergestes_atlanticus 0.0296403
## 63         Deosergestes_henseni            Sergia_tenuiremis 1.0000000
## 64    Gardinerosergia_splendens        Parasergestes_armatus 1.0000000
## 65    Gardinerosergia_splendens        Parasergestes_vigilax 0.0052508
## 66    Gardinerosergia_splendens         Phorcosergia_grandis 1.0000000
## 67    Gardinerosergia_splendens        Robustasergia_robusta 1.0000000
## 68    Gardinerosergia_splendens        Robustosergia_regalis 1.0000000
## 69    Gardinerosergia_splendens         Sergestes_atlanticus 0.1056249
## 70    Gardinerosergia_splendens            Sergia_tenuiremis 1.0000000
## 71        Parasergestes_armatus        Parasergestes_vigilax 0.0402557
## 72        Parasergestes_armatus         Phorcosergia_grandis 1.0000000
## 73        Parasergestes_armatus        Robustasergia_robusta 1.0000000
## 74        Parasergestes_armatus        Robustosergia_regalis 1.0000000
## 75        Parasergestes_armatus         Sergestes_atlanticus 0.5325877
## 76        Parasergestes_armatus            Sergia_tenuiremis 1.0000000
## 77        Parasergestes_vigilax         Phorcosergia_grandis 0.0023498
## 78        Parasergestes_vigilax        Robustasergia_robusta 0.0023706
## 79        Parasergestes_vigilax        Robustosergia_regalis 0.0016028
## 80        Parasergestes_vigilax         Sergestes_atlanticus 1.0000000
## 81        Parasergestes_vigilax            Sergia_tenuiremis 0.0025428
## 82         Phorcosergia_grandis        Robustasergia_robusta 1.0000000
## 83         Phorcosergia_grandis        Robustosergia_regalis 1.0000000
## 84         Phorcosergia_grandis         Sergestes_atlanticus 0.0579011
## 85         Phorcosergia_grandis            Sergia_tenuiremis 1.0000000
## 86        Robustasergia_robusta        Robustosergia_regalis 1.0000000
## 87        Robustasergia_robusta         Sergestes_atlanticus 0.0527121
## 88        Robustasergia_robusta            Sergia_tenuiremis 1.0000000
## 89        Robustosergia_regalis         Sergestes_atlanticus 0.0398191
## 90        Robustosergia_regalis            Sergia_tenuiremis 1.0000000
## 91         Sergestes_atlanticus            Sergia_tenuiremis 0.0504742
##        TestStat
## 1  1.962212e+00
## 2  9.226329e+00
## 3  1.043237e+01
## 4  1.571797e+01
## 5  1.472380e+01
## 6  1.245939e+01
## 7  8.972336e+00
## 8  3.803004e-01
## 9  1.342220e+01
## 10 1.372312e+01
## 11 1.411991e+01
## 12 2.222012e-01
## 13 1.389088e+01
## 14 3.508068e+00
## 15 4.594967e+00
## 16 1.077842e+01
## 17 9.396908e+00
## 18 6.573564e+00
## 19 3.374243e+00
## 20 3.924426e+00
## 21 8.390496e+00
## 22 8.028775e+00
## 23 9.130668e+00
## 24 8.872167e-01
## 25 8.010586e+00
## 26 3.390306e-03
## 27 2.657547e+00
## 28 1.879411e+00
## 29 8.991260e-01
## 30 2.114558e-03
## 31 1.266845e+01
## 32 7.471125e-01
## 33 1.399926e+00
## 34 1.263374e+00
## 35 7.202759e+00
## 36 1.971425e+00
## 37 3.517454e+00
## 38 2.443692e+00
## 39 1.050667e+00
## 40 1.587944e-06
## 41 1.420276e+01
## 42 1.176573e+00
## 43 1.729630e+00
## 44 1.864296e+00
## 45 8.439223e+00
## 46 2.256888e+00
## 47 5.038606e-02
## 48 2.215773e-01
## 49 2.213473e+00
## 50 2.009085e+01
## 51 1.558854e+00
## 52 1.172183e-01
## 53 6.377370e-01
## 54 1.396222e+01
## 55 2.728885e-02
## 56 7.176646e-02
## 57 1.560527e+00
## 58 1.892932e+01
## 59 8.193544e-01
## 60 1.578903e-02
## 61 2.516278e-01
## 62 1.288861e+01
## 63 1.086722e-01
## 64 7.497359e-01
## 65 1.617170e+01
## 66 1.761039e-01
## 67 2.145849e-02
## 68 1.364838e-02
## 69 1.045082e+01
## 70 2.726872e-01
## 71 1.230636e+01
## 72 5.417472e-01
## 73 1.166447e+00
## 74 9.775638e-01
## 75 6.962666e+00
## 76 1.736730e+00
## 77 1.770075e+01
## 78 1.768397e+01
## 79 1.842986e+01
## 80 1.161856e+00
## 81 1.755051e+01
## 82 4.600589e-01
## 83 2.486339e-01
## 84 1.161214e+01
## 85 1.016510e+00
## 86 1.038193e-01
## 87 1.179186e+01
## 88 1.758664e-01
## 89 1.232713e+01
## 90 5.282953e-01
## 91 1.187480e+01
## 
## ------------------------------------------------------------
## Coefficients by group in variable "genus_species"
## 
## Group: Allosergestes_pectinatus 
##             elevation    slope
## estimate    -2.908342 2.009753
## lower limit -3.816721 1.427844
## upper limit -1.999963 2.828814
## 
## H0 : variables uncorrelated.
## R-squared : 0.6929183 
## P-value : 0.00022059 
## 
## Group: Allosergestes_sargassi 
##             elevation    slope
## estimate    -2.212583 1.474440
## lower limit -2.806943 1.104000
## upper limit -1.618222 1.969178
## 
## H0 : variables uncorrelated.
## R-squared : 0.246605 
## P-value : 0.0015179 
## 
## Group: Challengerosergia_hansjacobi 
##             elevation     slope
## estimate    -1.516521 1.0399481
## lower limit -1.879267 0.8284462
## upper limit -1.153776 1.3054464
## 
## H0 : variables uncorrelated.
## R-squared : 0.4712721 
## P-value : 3.731e-07 
## 
## Group: Challengerosergia_talismani 
##             elevation     slope
## estimate    -1.564958 1.0316438
## lower limit -1.815606 0.8835032
## upper limit -1.314310 1.2046238
## 
## H0 : variables uncorrelated.
## R-squared : 0.825775 
## P-value : 6.5338e-13 
## 
## Group: Deosergestes_corniculum 
##             elevation     slope
## estimate    -1.337193 0.8105777
## lower limit -1.618827 0.6598417
## upper limit -1.055559 0.9957482
## 
## H0 : variables uncorrelated.
## R-squared : 0.8471468 
## P-value : 6.2959e-08 
## 
## Group: Deosergestes_henseni 
##             elevation     slope
## estimate    -1.327868 0.8376931
## lower limit -1.616552 0.6763970
## upper limit -1.039183 1.0374523
## 
## H0 : variables uncorrelated.
## R-squared : 0.3519041 
## P-value : 9.2412e-07 
## 
## Group: Gardinerosergia_splendens 
##              elevation     slope
## estimate    -1.2559422 0.8778027
## lower limit -1.6211316 0.6687602
## upper limit -0.8907528 1.1521882
## 
## H0 : variables uncorrelated.
## R-squared : 0.2574859 
## P-value : 0.00060377 
## 
## Group: Parasergestes_armatus 
##             elevation     slope
## estimate    -1.660526 1.0318360
## lower limit -2.051439 0.8006561
## upper limit -1.269613 1.3297664
## 
## H0 : variables uncorrelated.
## R-squared : 0.5613314 
## P-value : 1.9015e-06 
## 
## Group: Parasergestes_vigilax 
##             elevation    slope
## estimate    -3.256619 2.328849
## lower limit -4.339295 1.631234
## upper limit -2.173943 3.324807
## 
## H0 : variables uncorrelated.
## R-squared : 0.5634777 
## P-value : 0.00051661 
## 
## Group: Phorcosergia_grandis 
##             elevation     slope
## estimate    -1.463079 0.9333209
## lower limit -1.625863 0.8440017
## upper limit -1.300296 1.0320925
## 
## H0 : variables uncorrelated.
## R-squared : 0.8768029 
## P-value : < 2.22e-16 
## 
## Group: Robustasergia_robusta 
##              elevation     slope
## estimate    -1.2367827 0.8547738
## lower limit -1.6029678 0.6709735
## upper limit -0.8705977 1.0889227
## 
## H0 : variables uncorrelated.
## R-squared : 0.7097726 
## P-value : 4.5859e-07 
## 
## Group: Robustosergia_regalis 
##             elevation     slope
## estimate    -1.317150 0.8937637
## lower limit -1.526402 0.7755937
## upper limit -1.107898 1.0299381
## 
## H0 : variables uncorrelated.
## R-squared : 0.8227993 
## P-value : 4.2961e-15 
## 
## Group: Sergestes_atlanticus 
##             elevation    slope
## estimate    -2.629861 1.805583
## lower limit -3.474390 1.301400
## upper limit -1.785333 2.505095
## 
## H0 : variables uncorrelated.
## R-squared : 0.6905247 
## P-value : 0.00012413 
## 
## Group: Sergia_tenuiremis 
##              elevation     slope
## estimate    -1.1996991 0.7854599
## lower limit -1.6293146 0.5673032
## upper limit -0.7700836 1.0875090
## 
## H0 : variables uncorrelated.
## R-squared : 0.3252106 
## P-value : 0.0015327
#save coefficients of fits as object
cc_species <- data.frame(coef(sma_species))

This figure shows the developmental allometry of all the species with sufficient sampling (n = 14). It’s interactive, so you can click on different species in the legend to hide or show those points. That can be helpful for exploring patterns in a big mess of data. To me, it looks a bit like smaller-bodied species have higher slopes? That’s something we could test for if we wanted to. Here, it looks like Sergestes atlanticus, Allosergestes_pectinatus, and Parasergestes vigilax have different (higher) slopes compared to the other sergestids, and the other species differ a bit in slopes but more so in intercepts.

#make shape pallette
shapes.sp <- c("Allosergestes pectinatus" = 21,
              "Allosergestes sargassi" = 22, 
              "Challengerosergia hansjacobi" = 23, 
              "Challengerosergia talismani" = 24,
              "Deosergestes corniculum" = 25 , 
              "Deosergestes henseni" = 21, 
              "Gardinerosergia splendens" = 22, 
              "Parasergestes armatus" = 23, 
              "Parasergestes vigilax" = 24, 
              "Phorcosergia grandis" = 25, 
              "Robustasergia robusta" = 21 ,
              "Robustosergia regalis" = 22,
              "Sergestes atlanticus" = 23,
              "Sergia tenuiremis" = 24)

#make color pallette
cols.sp <- c("Allosergestes pectinatus" = "#a6cee3",
              "Allosergestes sargassi" = "#1f78b4", 
              "Challengerosergia hansjacobi" = "#b2df8a", 
              "Challengerosergia talismani" = "#33a02c",
              "Deosergestes corniculum" = "#fb9a99" , 
              "Deosergestes henseni" = "#e31a1c", 
              "Gardinerosergia splendens" = "#fdbf6f", 
              "Parasergestes armatus" = "#ff7f00", 
              "Parasergestes vigilax" = "#cab2d6", 
              "Phorcosergia grandis" = "#6a3d9a", 
              "Robustasergia robusta" ="gray31" ,
              "Robustosergia regalis" = "black",
              "Sergestes atlanticus" = "maroon1",
              "Sergia tenuiremis" = "orchid1")

#set s limits by species for plots
limits <- specimens_sma %>%
  group_by(genus_species) %>%
  summarise(max = max(Body_Length, na.rm = TRUE),
            min = min(Body_Length, na.rm = TRUE))


#make plot
plot_species <- ggplot(specimens_sma, aes(x = Body_Length, y = Eye_Diameter, color =  species_text, shape = species_text, fill = species_text)) + 
  geom_point(size = 2.5, alpha = 0.7) + 
  scale_shape_manual(values = shapes.sp, name = "Species") + 
  scale_color_manual(values = cols.sp, name = "Species") +
  scale_fill_manual(values = cols.sp, name = "Species") +
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(), legend.text = element_text(face = "italic")) +
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Allosergestes_pectinatus")], xend = limits$max[which(limits$genus_species == "Allosergestes_pectinatus")], y = 10^(log10(limits$min[which(limits$genus_species == "Allosergestes_pectinatus")])*cc_species["Allosergestes_pectinatus", "slope"] + cc_species["Allosergestes_pectinatus", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Allosergestes_pectinatus")])*cc_species["Allosergestes_pectinatus", "slope"] + cc_species["Allosergestes_pectinatus", "elevation"])), color = cols.sp["Allosergestes pectinatus"]) + #Rana temporaria adult scaling
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Allosergestes_sargassi")], xend = limits$max[which(limits$genus_species == "Allosergestes_sargassi")], y = 10^(log10(limits$min[which(limits$genus_species == "Allosergestes_sargassi")])*cc_species["Allosergestes_sargassi", "slope"] + cc_species["Allosergestes_sargassi", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Allosergestes_sargassi")])*cc_species["Allosergestes_sargassi", "slope"] + cc_species["Allosergestes_sargassi", "elevation"])), color = cols.sp["Allosergestes sargassi"]) + #Allosergestes_sargassi adult scaling
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Challengerosergia_hansjacobi")], xend = limits$max[which(limits$genus_species == "Challengerosergia_hansjacobi")], y = 10^(log10(limits$min[which(limits$genus_species == "Challengerosergia_hansjacobi")])*cc_species["Challengerosergia_hansjacobi", "slope"] + cc_species["Challengerosergia_hansjacobi", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Challengerosergia_hansjacobi")])*cc_species["Challengerosergia_hansjacobi", "slope"] + cc_species["Challengerosergia_hansjacobi", "elevation"])), color = cols.sp["Challengerosergia hansjacobi"]) + #Challengerosergia_hansjacobi adult scaling
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Challengerosergia_talismani")], xend = limits$max[which(limits$genus_species == "Challengerosergia_talismani")], y = 10^(log10(limits$min[which(limits$genus_species == "Challengerosergia_talismani")])*cc_species["Challengerosergia_talismani", "slope"] + cc_species["Challengerosergia_talismani", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Challengerosergia_talismani")])*cc_species["Challengerosergia_talismani", "slope"] + cc_species["Challengerosergia_talismani", "elevation"])), color = cols.sp["Challengerosergia talismani"]) + #Challengerosergia_talismani adult scaling
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Deosergestes_corniculum")], xend = limits$max[which(limits$genus_species == "Deosergestes_corniculum")], y = 10^(log10(limits$min[which(limits$genus_species == "Deosergestes_corniculum")])*cc_species["Deosergestes_corniculum", "slope"] + cc_species["Deosergestes_corniculum", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Deosergestes_corniculum")])*cc_species["Deosergestes_corniculum", "slope"] + cc_species["Deosergestes_corniculum", "elevation"])), color = cols.sp["Deosergestes corniculum"]) + #Deosergestes_corniculum adult scaling
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Deosergestes_henseni")], xend = limits$max[which(limits$genus_species == "Deosergestes_henseni")], y = 10^(log10(limits$min[which(limits$genus_species == "Deosergestes_henseni")])*cc_species["Deosergestes_henseni", "slope"] + cc_species["Deosergestes_henseni", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Deosergestes_henseni")])*cc_species["Deosergestes_henseni", "slope"] + cc_species["Deosergestes_henseni", "elevation"])), color = cols.sp["Deosergestes henseni"]) + #Deosergestes henseni adult scaling
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Gardinerosergia_splendens")], xend = limits$max[which(limits$genus_species == "Gardinerosergia_splendens")], y = 10^(log10(limits$min[which(limits$genus_species == "Gardinerosergia_splendens")])*cc_species["Gardinerosergia_splendens", "slope"] + cc_species["Gardinerosergia_splendens", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Gardinerosergia_splendens")])*cc_species["Gardinerosergia_splendens", "slope"] + cc_species["Gardinerosergia_splendens", "elevation"])), color = cols.sp["Gardinerosergia splendens"]) + #Gardinerosergia splendens adult scaling
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Parasergestes_armatus")], xend = limits$max[which(limits$genus_species == "Parasergestes_armatus")], y = 10^(log10(limits$min[which(limits$genus_species == "Parasergestes_armatus")])*cc_species["Parasergestes_armatus", "slope"] + cc_species["Parasergestes_armatus", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Parasergestes_armatus")])*cc_species["Parasergestes_armatus", "slope"] + cc_species["Parasergestes_armatus", "elevation"])), color = cols.sp["Parasergestes armatus"]) + #Parasergestes armatus adult scaling
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Parasergestes_vigilax")], xend = limits$max[which(limits$genus_species == "Parasergestes_vigilax")], y = 10^(log10(limits$min[which(limits$genus_species == "Parasergestes_vigilax")])*cc_species["Parasergestes_vigilax", "slope"] + cc_species["Parasergestes_vigilax", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Parasergestes_vigilax")])*cc_species["Parasergestes_vigilax", "slope"] + cc_species["Parasergestes_vigilax", "elevation"])), color = cols.sp["Parasergestes vigilax"]) + #Parasergestes vigilax adult scaling
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Phorcosergia_grandis")], xend = limits$max[which(limits$genus_species == "Phorcosergia_grandis")], y = 10^(log10(limits$min[which(limits$genus_species == "Phorcosergia_grandis")])*cc_species["Phorcosergia_grandis", "slope"] + cc_species["Phorcosergia_grandis", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Phorcosergia_grandis")])*cc_species["Phorcosergia_grandis", "slope"] + cc_species["Phorcosergia_grandis", "elevation"])), color = cols.sp["Hyla meridiondalis"]) + #Phorcosergia grandis adult scaling
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Robustasergia_robusta")], xend = limits$max[which(limits$genus_species == "Robustasergia_robusta")], y = 10^(log10(limits$min[which(limits$genus_species == "Robustasergia_robusta")])*cc_species["Robustasergia_robusta", "slope"] + cc_species["Robustasergia_robusta", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Robustasergia_robusta")])*cc_species["Robustasergia_robusta", "slope"] + cc_species["Robustasergia_robusta", "elevation"])), color = cols.sp["Robustasergia robusta"]) + #Robustasergia robusta adult scaling
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Robustosergia_regalis")], xend = limits$max[which(limits$genus_species == "Robustosergia_regalis")], y = 10^(log10(limits$min[which(limits$genus_species == "Robustosergia_regalis")])*cc_species["Robustosergia_regalis", "slope"] + cc_species["Robustosergia_regalis", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Robustosergia_regalis")])*cc_species["Robustosergia_regalis", "slope"] + cc_species["Robustosergia_regalis", "elevation"])), color = cols.sp["Robustosergia regalis"]) +
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Sergestes_atlanticus")], xend = limits$max[which(limits$genus_species == "Sergestes_atlanticus")], y = 10^(log10(limits$min[which(limits$genus_species == "Sergestes_atlanticus")])*cc_species["Sergestes_atlanticus", "slope"] + cc_species["Sergestes_atlanticus", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Sergestes_atlanticus")])*cc_species["Sergestes_atlanticus", "slope"] + cc_species["Sergestes_atlanticus", "elevation"])), color = cols.sp["Sergestes atlanticus"]) +
  geom_segment(aes(x = limits$min[which(limits$genus_species == "Sergia_tenuiremis")], xend = limits$max[which(limits$genus_species == "Sergia_tenuiremis")], y = 10^(log10(limits$min[which(limits$genus_species == "Sergia_tenuiremis")])*cc_species["Sergia_tenuiremis", "slope"] + cc_species["Sergia_tenuiremis", "elevation"]), yend = 10^(log10(limits$max[which(limits$genus_species == "Sergia_tenuiremis")])*cc_species["Sergia_tenuiremis", "slope"] + cc_species["Sergia_tenuiremis", "elevation"])), color = cols.sp["Sergia tenuiremis"])


#interactive plot
ggplotly(plot_species)
#export figure
png("../Figures/dev-species.png", width = 12, height = 12, units = "in", res = 500)
plot_species
dev.off()

Test for effect of preservation

As one additional possible test to see if preservation type affects allometric fits, I’ve included a linear model here that has all of the specimen data and fits a regression of eye size vs. body size with both species and preservation type as covariates (and looks for an interaction effect). Specifically, this tests whether preservation type (fresh, ethanol, paraformaldehyde) has an effect on eye size, and whether this differs across species.

# Model for developmental scaling across specimens * species * preservation
lm_pres <- lm(formula = Eye_Diameter ~ Body_Length * genus_species * Preservation, data = specimens_sma) 

#view regression estimates and t-tests
summary(lm_pres)
## 
## Call:
## lm(formula = Eye_Diameter ~ Body_Length * genus_species * Preservation, 
##     data = specimens_sma)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.45611 -0.07657 -0.00360  0.06465  0.84951 
## 
## Coefficients: (18 not defined because of singularities)
##                                                                                      Estimate
## (Intercept)                                                                        -3.915e-02
## Body_Length                                                                         2.373e-02
## genus_speciesAllosergestes_sargassi                                                -9.075e-02
## genus_speciesChallengerosergia_hansjacobi                                           3.981e-01
## genus_speciesChallengerosergia_talismani                                            1.361e-01
## genus_speciesDeosergestes_corniculum                                                3.035e-01
## genus_speciesDeosergestes_henseni                                                   4.553e-01
## genus_speciesGardinerosergia_splendens                                              4.834e-01
## genus_speciesParasergestes_armatus                                                  1.287e-01
## genus_speciesParasergestes_vigilax                                                  1.266e-01
## genus_speciesPhorcosergia_grandis                                                   2.190e-01
## genus_speciesRobustasergia_robusta                                                  5.868e-01
## genus_speciesRobustosergia_regalis                                                  3.460e-01
## genus_speciesSergestes_atlanticus                                                   1.142e-01
## genus_speciesSergia_tenuiremis                                                      1.738e+00
## Preservationfresh                                                                  -1.001e+00
## Preservationparaformaldehyde                                                       -3.007e-01
## Body_Length:genus_speciesAllosergestes_sargassi                                     4.758e-03
## Body_Length:genus_speciesChallengerosergia_hansjacobi                               1.057e-03
## Body_Length:genus_speciesChallengerosergia_talismani                                4.224e-03
## Body_Length:genus_speciesDeosergestes_corniculum                                   -7.388e-03
## Body_Length:genus_speciesDeosergestes_henseni                                      -1.002e-02
## Body_Length:genus_speciesGardinerosergia_splendens                                 -6.092e-03
## Body_Length:genus_speciesParasergestes_armatus                                     -5.120e-03
## Body_Length:genus_speciesParasergestes_vigilax                                      1.053e-03
## Body_Length:genus_speciesPhorcosergia_grandis                                      -9.112e-04
## Body_Length:genus_speciesRobustasergia_robusta                                     -1.162e-03
## Body_Length:genus_speciesRobustosergia_regalis                                      1.264e-03
## Body_Length:genus_speciesSergestes_atlanticus                                       3.619e-05
## Body_Length:genus_speciesSergia_tenuiremis                                         -6.029e-02
## Body_Length:Preservationfresh                                                       4.895e-02
## Body_Length:Preservationparaformaldehyde                                            2.338e-02
## genus_speciesAllosergestes_sargassi:Preservationfresh                               1.595e+00
## genus_speciesChallengerosergia_hansjacobi:Preservationfresh                         1.146e+00
## genus_speciesChallengerosergia_talismani:Preservationfresh                          1.132e+00
## genus_speciesDeosergestes_corniculum:Preservationfresh                              1.329e+00
## genus_speciesDeosergestes_henseni:Preservationfresh                                 1.155e+00
## genus_speciesGardinerosergia_splendens:Preservationfresh                            1.320e+00
## genus_speciesParasergestes_armatus:Preservationfresh                                1.433e+00
## genus_speciesParasergestes_vigilax:Preservationfresh                               -1.019e-01
## genus_speciesPhorcosergia_grandis:Preservationfresh                                 7.928e-01
## genus_speciesRobustasergia_robusta:Preservationfresh                                       NA
## genus_speciesRobustosergia_regalis:Preservationfresh                                1.299e+00
## genus_speciesSergestes_atlanticus:Preservationfresh                                 1.020e+00
## genus_speciesSergia_tenuiremis:Preservationfresh                                           NA
## genus_speciesAllosergestes_sargassi:Preservationparaformaldehyde                    2.470e+00
## genus_speciesChallengerosergia_hansjacobi:Preservationparaformaldehyde              5.159e-01
## genus_speciesChallengerosergia_talismani:Preservationparaformaldehyde              -4.247e-01
## genus_speciesDeosergestes_corniculum:Preservationparaformaldehyde                          NA
## genus_speciesDeosergestes_henseni:Preservationparaformaldehyde                             NA
## genus_speciesGardinerosergia_splendens:Preservationparaformaldehyde                        NA
## genus_speciesParasergestes_armatus:Preservationparaformaldehyde                     2.775e-01
## genus_speciesParasergestes_vigilax:Preservationparaformaldehyde                    -6.141e-01
## genus_speciesPhorcosergia_grandis:Preservationparaformaldehyde                      2.362e-01
## genus_speciesRobustasergia_robusta:Preservationparaformaldehyde                            NA
## genus_speciesRobustosergia_regalis:Preservationparaformaldehyde                            NA
## genus_speciesSergestes_atlanticus:Preservationparaformaldehyde                      7.713e-01
## genus_speciesSergia_tenuiremis:Preservationparaformaldehyde                                NA
## Body_Length:genus_speciesAllosergestes_sargassi:Preservationfresh                  -6.751e-02
## Body_Length:genus_speciesChallengerosergia_hansjacobi:Preservationfresh            -5.611e-02
## Body_Length:genus_speciesChallengerosergia_talismani:Preservationfresh             -5.357e-02
## Body_Length:genus_speciesDeosergestes_corniculum:Preservationfresh                 -5.439e-02
## Body_Length:genus_speciesDeosergestes_henseni:Preservationfresh                    -5.061e-02
## Body_Length:genus_speciesGardinerosergia_splendens:Preservationfresh               -5.263e-02
## Body_Length:genus_speciesParasergestes_armatus:Preservationfresh                   -5.853e-02
## Body_Length:genus_speciesParasergestes_vigilax:Preservationfresh                           NA
## Body_Length:genus_speciesPhorcosergia_grandis:Preservationfresh                    -4.444e-02
## Body_Length:genus_speciesRobustasergia_robusta:Preservationfresh                           NA
## Body_Length:genus_speciesRobustosergia_regalis:Preservationfresh                   -5.457e-02
## Body_Length:genus_speciesSergestes_atlanticus:Preservationfresh                    -4.433e-02
## Body_Length:genus_speciesSergia_tenuiremis:Preservationfresh                               NA
## Body_Length:genus_speciesAllosergestes_sargassi:Preservationparaformaldehyde       -1.084e-01
## Body_Length:genus_speciesChallengerosergia_hansjacobi:Preservationparaformaldehyde -3.198e-02
## Body_Length:genus_speciesChallengerosergia_talismani:Preservationparaformaldehyde          NA
## Body_Length:genus_speciesDeosergestes_corniculum:Preservationparaformaldehyde              NA
## Body_Length:genus_speciesDeosergestes_henseni:Preservationparaformaldehyde                 NA
## Body_Length:genus_speciesGardinerosergia_splendens:Preservationparaformaldehyde            NA
## Body_Length:genus_speciesParasergestes_armatus:Preservationparaformaldehyde        -1.455e-02
## Body_Length:genus_speciesParasergestes_vigilax:Preservationparaformaldehyde         2.378e-02
## Body_Length:genus_speciesPhorcosergia_grandis:Preservationparaformaldehyde         -1.958e-02
## Body_Length:genus_speciesRobustasergia_robusta:Preservationparaformaldehyde                NA
## Body_Length:genus_speciesRobustosergia_regalis:Preservationparaformaldehyde                NA
## Body_Length:genus_speciesSergestes_atlanticus:Preservationparaformaldehyde         -2.915e-02
## Body_Length:genus_speciesSergia_tenuiremis:Preservationparaformaldehyde                    NA
##                                                                                    Std. Error
## (Intercept)                                                                         4.056e-01
## Body_Length                                                                         2.193e-02
## genus_speciesAllosergestes_sargassi                                                 6.612e-01
## genus_speciesChallengerosergia_hansjacobi                                           4.765e-01
## genus_speciesChallengerosergia_talismani                                            4.291e-01
## genus_speciesDeosergestes_corniculum                                                4.300e-01
## genus_speciesDeosergestes_henseni                                                   4.709e-01
## genus_speciesGardinerosergia_splendens                                              4.244e-01
## genus_speciesParasergestes_armatus                                                  4.334e-01
## genus_speciesParasergestes_vigilax                                                  5.193e-01
## genus_speciesPhorcosergia_grandis                                                   4.145e-01
## genus_speciesRobustasergia_robusta                                                  4.288e-01
## genus_speciesRobustosergia_regalis                                                  4.202e-01
## genus_speciesSergestes_atlanticus                                                   7.168e-01
## genus_speciesSergia_tenuiremis                                                      9.428e-01
## Preservationfresh                                                                   1.016e+00
## Preservationparaformaldehyde                                                        5.604e-01
## Body_Length:genus_speciesAllosergestes_sargassi                                     3.341e-02
## Body_Length:genus_speciesChallengerosergia_hansjacobi                               2.311e-02
## Body_Length:genus_speciesChallengerosergia_talismani                                2.221e-02
## Body_Length:genus_speciesDeosergestes_corniculum                                    2.211e-02
## Body_Length:genus_speciesDeosergestes_henseni                                       2.282e-02
## Body_Length:genus_speciesGardinerosergia_splendens                                  2.211e-02
## Body_Length:genus_speciesParasergestes_armatus                                      2.255e-02
## Body_Length:genus_speciesParasergestes_vigilax                                      2.768e-02
## Body_Length:genus_speciesPhorcosergia_grandis                                       2.197e-02
## Body_Length:genus_speciesRobustasergia_robusta                                      2.206e-02
## Body_Length:genus_speciesRobustosergia_regalis                                      2.206e-02
## Body_Length:genus_speciesSergestes_atlanticus                                       3.379e-02
## Body_Length:genus_speciesSergia_tenuiremis                                          4.341e-02
## Body_Length:Preservationfresh                                                       4.854e-02
## Body_Length:Preservationparaformaldehyde                                            2.867e-02
## genus_speciesAllosergestes_sargassi:Preservationfresh                               1.160e+00
## genus_speciesChallengerosergia_hansjacobi:Preservationfresh                         3.296e+00
## genus_speciesChallengerosergia_talismani:Preservationfresh                          1.051e+00
## genus_speciesDeosergestes_corniculum:Preservationfresh                              2.860e+00
## genus_speciesDeosergestes_henseni:Preservationfresh                                 1.049e+00
## genus_speciesGardinerosergia_splendens:Preservationfresh                            1.046e+00
## genus_speciesParasergestes_armatus:Preservationfresh                                1.139e+00
## genus_speciesParasergestes_vigilax:Preservationfresh                                2.890e-01
## genus_speciesPhorcosergia_grandis:Preservationfresh                                 1.037e+00
## genus_speciesRobustasergia_robusta:Preservationfresh                                       NA
## genus_speciesRobustosergia_regalis:Preservationfresh                                1.113e+00
## genus_speciesSergestes_atlanticus:Preservationfresh                                 1.336e+00
## genus_speciesSergia_tenuiremis:Preservationfresh                                           NA
## genus_speciesAllosergestes_sargassi:Preservationparaformaldehyde                    1.653e+00
## genus_speciesChallengerosergia_hansjacobi:Preservationparaformaldehyde              7.155e-01
## genus_speciesChallengerosergia_talismani:Preservationparaformaldehyde               6.984e-01
## genus_speciesDeosergestes_corniculum:Preservationparaformaldehyde                          NA
## genus_speciesDeosergestes_henseni:Preservationparaformaldehyde                             NA
## genus_speciesGardinerosergia_splendens:Preservationparaformaldehyde                        NA
## genus_speciesParasergestes_armatus:Preservationparaformaldehyde                     6.464e-01
## genus_speciesParasergestes_vigilax:Preservationparaformaldehyde                     7.692e-01
## genus_speciesPhorcosergia_grandis:Preservationparaformaldehyde                      7.375e-01
## genus_speciesRobustasergia_robusta:Preservationparaformaldehyde                            NA
## genus_speciesRobustosergia_regalis:Preservationparaformaldehyde                            NA
## genus_speciesSergestes_atlanticus:Preservationparaformaldehyde                      1.136e+00
## genus_speciesSergia_tenuiremis:Preservationparaformaldehyde                                NA
## Body_Length:genus_speciesAllosergestes_sargassi:Preservationfresh                   5.526e-02
## Body_Length:genus_speciesChallengerosergia_hansjacobi:Preservationfresh             1.032e-01
## Body_Length:genus_speciesChallengerosergia_talismani:Preservationfresh              4.911e-02
## Body_Length:genus_speciesDeosergestes_corniculum:Preservationfresh                  7.427e-02
## Body_Length:genus_speciesDeosergestes_henseni:Preservationfresh                     4.900e-02
## Body_Length:genus_speciesGardinerosergia_splendens:Preservationfresh                4.912e-02
## Body_Length:genus_speciesParasergestes_armatus:Preservationfresh                    5.078e-02
## Body_Length:genus_speciesParasergestes_vigilax:Preservationfresh                           NA
## Body_Length:genus_speciesPhorcosergia_grandis:Preservationfresh                     4.873e-02
## Body_Length:genus_speciesRobustasergia_robusta:Preservationfresh                           NA
## Body_Length:genus_speciesRobustosergia_regalis:Preservationfresh                    4.951e-02
## Body_Length:genus_speciesSergestes_atlanticus:Preservationfresh                     6.004e-02
## Body_Length:genus_speciesSergia_tenuiremis:Preservationfresh                               NA
## Body_Length:genus_speciesAllosergestes_sargassi:Preservationparaformaldehyde        7.551e-02
## Body_Length:genus_speciesChallengerosergia_hansjacobi:Preservationparaformaldehyde  3.186e-02
## Body_Length:genus_speciesChallengerosergia_talismani:Preservationparaformaldehyde          NA
## Body_Length:genus_speciesDeosergestes_corniculum:Preservationparaformaldehyde              NA
## Body_Length:genus_speciesDeosergestes_henseni:Preservationparaformaldehyde                 NA
## Body_Length:genus_speciesGardinerosergia_splendens:Preservationparaformaldehyde            NA
## Body_Length:genus_speciesParasergestes_armatus:Preservationparaformaldehyde         3.088e-02
## Body_Length:genus_speciesParasergestes_vigilax:Preservationparaformaldehyde         3.917e-02
## Body_Length:genus_speciesPhorcosergia_grandis:Preservationparaformaldehyde          3.093e-02
## Body_Length:genus_speciesRobustasergia_robusta:Preservationparaformaldehyde                NA
## Body_Length:genus_speciesRobustosergia_regalis:Preservationparaformaldehyde                NA
## Body_Length:genus_speciesSergestes_atlanticus:Preservationparaformaldehyde          4.786e-02
## Body_Length:genus_speciesSergia_tenuiremis:Preservationparaformaldehyde                    NA
##                                                                                    t value
## (Intercept)                                                                         -0.097
## Body_Length                                                                          1.082
## genus_speciesAllosergestes_sargassi                                                 -0.137
## genus_speciesChallengerosergia_hansjacobi                                            0.835
## genus_speciesChallengerosergia_talismani                                             0.317
## genus_speciesDeosergestes_corniculum                                                 0.706
## genus_speciesDeosergestes_henseni                                                    0.967
## genus_speciesGardinerosergia_splendens                                               1.139
## genus_speciesParasergestes_armatus                                                   0.297
## genus_speciesParasergestes_vigilax                                                   0.244
## genus_speciesPhorcosergia_grandis                                                    0.528
## genus_speciesRobustasergia_robusta                                                   1.369
## genus_speciesRobustosergia_regalis                                                   0.823
## genus_speciesSergestes_atlanticus                                                    0.159
## genus_speciesSergia_tenuiremis                                                       1.843
## Preservationfresh                                                                   -0.986
## Preservationparaformaldehyde                                                        -0.537
## Body_Length:genus_speciesAllosergestes_sargassi                                      0.142
## Body_Length:genus_speciesChallengerosergia_hansjacobi                                0.046
## Body_Length:genus_speciesChallengerosergia_talismani                                 0.190
## Body_Length:genus_speciesDeosergestes_corniculum                                    -0.334
## Body_Length:genus_speciesDeosergestes_henseni                                       -0.439
## Body_Length:genus_speciesGardinerosergia_splendens                                  -0.275
## Body_Length:genus_speciesParasergestes_armatus                                      -0.227
## Body_Length:genus_speciesParasergestes_vigilax                                       0.038
## Body_Length:genus_speciesPhorcosergia_grandis                                       -0.041
## Body_Length:genus_speciesRobustasergia_robusta                                      -0.053
## Body_Length:genus_speciesRobustosergia_regalis                                       0.057
## Body_Length:genus_speciesSergestes_atlanticus                                        0.001
## Body_Length:genus_speciesSergia_tenuiremis                                          -1.389
## Body_Length:Preservationfresh                                                        1.008
## Body_Length:Preservationparaformaldehyde                                             0.815
## genus_speciesAllosergestes_sargassi:Preservationfresh                                1.376
## genus_speciesChallengerosergia_hansjacobi:Preservationfresh                          0.348
## genus_speciesChallengerosergia_talismani:Preservationfresh                           1.078
## genus_speciesDeosergestes_corniculum:Preservationfresh                               0.465
## genus_speciesDeosergestes_henseni:Preservationfresh                                  1.102
## genus_speciesGardinerosergia_splendens:Preservationfresh                             1.261
## genus_speciesParasergestes_armatus:Preservationfresh                                 1.258
## genus_speciesParasergestes_vigilax:Preservationfresh                                -0.353
## genus_speciesPhorcosergia_grandis:Preservationfresh                                  0.765
## genus_speciesRobustasergia_robusta:Preservationfresh                                    NA
## genus_speciesRobustosergia_regalis:Preservationfresh                                 1.168
## genus_speciesSergestes_atlanticus:Preservationfresh                                  0.764
## genus_speciesSergia_tenuiremis:Preservationfresh                                        NA
## genus_speciesAllosergestes_sargassi:Preservationparaformaldehyde                     1.494
## genus_speciesChallengerosergia_hansjacobi:Preservationparaformaldehyde               0.721
## genus_speciesChallengerosergia_talismani:Preservationparaformaldehyde               -0.608
## genus_speciesDeosergestes_corniculum:Preservationparaformaldehyde                       NA
## genus_speciesDeosergestes_henseni:Preservationparaformaldehyde                          NA
## genus_speciesGardinerosergia_splendens:Preservationparaformaldehyde                     NA
## genus_speciesParasergestes_armatus:Preservationparaformaldehyde                      0.429
## genus_speciesParasergestes_vigilax:Preservationparaformaldehyde                     -0.798
## genus_speciesPhorcosergia_grandis:Preservationparaformaldehyde                       0.320
## genus_speciesRobustasergia_robusta:Preservationparaformaldehyde                         NA
## genus_speciesRobustosergia_regalis:Preservationparaformaldehyde                         NA
## genus_speciesSergestes_atlanticus:Preservationparaformaldehyde                       0.679
## genus_speciesSergia_tenuiremis:Preservationparaformaldehyde                             NA
## Body_Length:genus_speciesAllosergestes_sargassi:Preservationfresh                   -1.222
## Body_Length:genus_speciesChallengerosergia_hansjacobi:Preservationfresh             -0.544
## Body_Length:genus_speciesChallengerosergia_talismani:Preservationfresh              -1.091
## Body_Length:genus_speciesDeosergestes_corniculum:Preservationfresh                  -0.732
## Body_Length:genus_speciesDeosergestes_henseni:Preservationfresh                     -1.033
## Body_Length:genus_speciesGardinerosergia_splendens:Preservationfresh                -1.072
## Body_Length:genus_speciesParasergestes_armatus:Preservationfresh                    -1.153
## Body_Length:genus_speciesParasergestes_vigilax:Preservationfresh                        NA
## Body_Length:genus_speciesPhorcosergia_grandis:Preservationfresh                     -0.912
## Body_Length:genus_speciesRobustasergia_robusta:Preservationfresh                        NA
## Body_Length:genus_speciesRobustosergia_regalis:Preservationfresh                    -1.102
## Body_Length:genus_speciesSergestes_atlanticus:Preservationfresh                     -0.738
## Body_Length:genus_speciesSergia_tenuiremis:Preservationfresh                            NA
## Body_Length:genus_speciesAllosergestes_sargassi:Preservationparaformaldehyde        -1.436
## Body_Length:genus_speciesChallengerosergia_hansjacobi:Preservationparaformaldehyde  -1.004
## Body_Length:genus_speciesChallengerosergia_talismani:Preservationparaformaldehyde       NA
## Body_Length:genus_speciesDeosergestes_corniculum:Preservationparaformaldehyde           NA
## Body_Length:genus_speciesDeosergestes_henseni:Preservationparaformaldehyde              NA
## Body_Length:genus_speciesGardinerosergia_splendens:Preservationparaformaldehyde         NA
## Body_Length:genus_speciesParasergestes_armatus:Preservationparaformaldehyde         -0.471
## Body_Length:genus_speciesParasergestes_vigilax:Preservationparaformaldehyde          0.607
## Body_Length:genus_speciesPhorcosergia_grandis:Preservationparaformaldehyde          -0.633
## Body_Length:genus_speciesRobustasergia_robusta:Preservationparaformaldehyde             NA
## Body_Length:genus_speciesRobustosergia_regalis:Preservationparaformaldehyde             NA
## Body_Length:genus_speciesSergestes_atlanticus:Preservationparaformaldehyde          -0.609
## Body_Length:genus_speciesSergia_tenuiremis:Preservationparaformaldehyde                 NA
##                                                                                    Pr(>|t|)
## (Intercept)                                                                          0.9232
## Body_Length                                                                          0.2800
## genus_speciesAllosergestes_sargassi                                                  0.8909
## genus_speciesChallengerosergia_hansjacobi                                            0.4040
## genus_speciesChallengerosergia_talismani                                             0.7513
## genus_speciesDeosergestes_corniculum                                                 0.4807
## genus_speciesDeosergestes_henseni                                                    0.3342
## genus_speciesGardinerosergia_splendens                                               0.2555
## genus_speciesParasergestes_armatus                                                   0.7666
## genus_speciesParasergestes_vigilax                                                   0.8075
## genus_speciesPhorcosergia_grandis                                                    0.5975
## genus_speciesRobustasergia_robusta                                                   0.1719
## genus_speciesRobustosergia_regalis                                                   0.4107
## genus_speciesSergestes_atlanticus                                                    0.8735
## genus_speciesSergia_tenuiremis                                                       0.0661
## Preservationfresh                                                                    0.3249
## Preservationparaformaldehyde                                                         0.5918
## Body_Length:genus_speciesAllosergestes_sargassi                                      0.8868
## Body_Length:genus_speciesChallengerosergia_hansjacobi                                0.9635
## Body_Length:genus_speciesChallengerosergia_talismani                                 0.8493
## Body_Length:genus_speciesDeosergestes_corniculum                                     0.7385
## Body_Length:genus_speciesDeosergestes_henseni                                        0.6608
## Body_Length:genus_speciesGardinerosergia_splendens                                   0.7831
## Body_Length:genus_speciesParasergestes_armatus                                       0.8205
## Body_Length:genus_speciesParasergestes_vigilax                                       0.9697
## Body_Length:genus_speciesPhorcosergia_grandis                                        0.9669
## Body_Length:genus_speciesRobustasergia_robusta                                       0.9580
## Body_Length:genus_speciesRobustosergia_regalis                                       0.9543
## Body_Length:genus_speciesSergestes_atlanticus                                        0.9991
## Body_Length:genus_speciesSergia_tenuiremis                                           0.1657
## Body_Length:Preservationfresh                                                        0.3139
## Body_Length:Preservationparaformaldehyde                                             0.4153
## genus_speciesAllosergestes_sargassi:Preservationfresh                                0.1697
## genus_speciesChallengerosergia_hansjacobi:Preservationfresh                          0.7283
## genus_speciesChallengerosergia_talismani:Preservationfresh                           0.2819
## genus_speciesDeosergestes_corniculum:Preservationfresh                               0.6423
## genus_speciesDeosergestes_henseni:Preservationfresh                                  0.2713
## genus_speciesGardinerosergia_splendens:Preservationfresh                             0.2081
## genus_speciesParasergestes_armatus:Preservationfresh                                 0.2092
## genus_speciesParasergestes_vigilax:Preservationfresh                                 0.7247
## genus_speciesPhorcosergia_grandis:Preservationfresh                                  0.4449
## genus_speciesRobustasergia_robusta:Preservationfresh                                     NA
## genus_speciesRobustosergia_regalis:Preservationfresh                                 0.2437
## genus_speciesSergestes_atlanticus:Preservationfresh                                  0.4455
## genus_speciesSergia_tenuiremis:Preservationfresh                                         NA
## genus_speciesAllosergestes_sargassi:Preservationparaformaldehyde                     0.1359
## genus_speciesChallengerosergia_hansjacobi:Preservationparaformaldehyde               0.4713
## genus_speciesChallengerosergia_talismani:Preservationparaformaldehyde                0.5435
## genus_speciesDeosergestes_corniculum:Preservationparaformaldehyde                        NA
## genus_speciesDeosergestes_henseni:Preservationparaformaldehyde                           NA
## genus_speciesGardinerosergia_splendens:Preservationparaformaldehyde                      NA
## genus_speciesParasergestes_armatus:Preservationparaformaldehyde                      0.6679
## genus_speciesParasergestes_vigilax:Preservationparaformaldehyde                      0.4251
## genus_speciesPhorcosergia_grandis:Preservationparaformaldehyde                       0.7490
## genus_speciesRobustasergia_robusta:Preservationparaformaldehyde                          NA
## genus_speciesRobustosergia_regalis:Preservationparaformaldehyde                          NA
## genus_speciesSergestes_atlanticus:Preservationparaformaldehyde                       0.4977
## genus_speciesSergia_tenuiremis:Preservationparaformaldehyde                              NA
## Body_Length:genus_speciesAllosergestes_sargassi:Preservationfresh                    0.2226
## Body_Length:genus_speciesChallengerosergia_hansjacobi:Preservationfresh              0.5869
## Body_Length:genus_speciesChallengerosergia_talismani:Preservationfresh               0.2760
## Body_Length:genus_speciesDeosergestes_corniculum:Preservationfresh                   0.4644
## Body_Length:genus_speciesDeosergestes_henseni:Preservationfresh                      0.3023
## Body_Length:genus_speciesGardinerosergia_splendens:Preservationfresh                 0.2846
## Body_Length:genus_speciesParasergestes_armatus:Preservationfresh                     0.2498
## Body_Length:genus_speciesParasergestes_vigilax:Preservationfresh                         NA
## Body_Length:genus_speciesPhorcosergia_grandis:Preservationfresh                      0.3623
## Body_Length:genus_speciesRobustasergia_robusta:Preservationfresh                         NA
## Body_Length:genus_speciesRobustosergia_regalis:Preservationfresh                     0.2711
## Body_Length:genus_speciesSergestes_atlanticus:Preservationfresh                      0.4608
## Body_Length:genus_speciesSergia_tenuiremis:Preservationfresh                             NA
## Body_Length:genus_speciesAllosergestes_sargassi:Preservationparaformaldehyde         0.1518
## Body_Length:genus_speciesChallengerosergia_hansjacobi:Preservationparaformaldehyde   0.3162
## Body_Length:genus_speciesChallengerosergia_talismani:Preservationparaformaldehyde        NA
## Body_Length:genus_speciesDeosergestes_corniculum:Preservationparaformaldehyde            NA
## Body_Length:genus_speciesDeosergestes_henseni:Preservationparaformaldehyde               NA
## Body_Length:genus_speciesGardinerosergia_splendens:Preservationparaformaldehyde          NA
## Body_Length:genus_speciesParasergestes_armatus:Preservationparaformaldehyde          0.6378
## Body_Length:genus_speciesParasergestes_vigilax:Preservationparaformaldehyde          0.5442
## Body_Length:genus_speciesPhorcosergia_grandis:Preservationparaformaldehyde           0.5272
## Body_Length:genus_speciesRobustasergia_robusta:Preservationparaformaldehyde              NA
## Body_Length:genus_speciesRobustosergia_regalis:Preservationparaformaldehyde              NA
## Body_Length:genus_speciesSergestes_atlanticus:Preservationparaformaldehyde           0.5429
## Body_Length:genus_speciesSergia_tenuiremis:Preservationparaformaldehyde                  NA
##                                                                                     
## (Intercept)                                                                         
## Body_Length                                                                         
## genus_speciesAllosergestes_sargassi                                                 
## genus_speciesChallengerosergia_hansjacobi                                           
## genus_speciesChallengerosergia_talismani                                            
## genus_speciesDeosergestes_corniculum                                                
## genus_speciesDeosergestes_henseni                                                   
## genus_speciesGardinerosergia_splendens                                              
## genus_speciesParasergestes_armatus                                                  
## genus_speciesParasergestes_vigilax                                                  
## genus_speciesPhorcosergia_grandis                                                   
## genus_speciesRobustasergia_robusta                                                  
## genus_speciesRobustosergia_regalis                                                  
## genus_speciesSergestes_atlanticus                                                   
## genus_speciesSergia_tenuiremis                                                     .
## Preservationfresh                                                                   
## Preservationparaformaldehyde                                                        
## Body_Length:genus_speciesAllosergestes_sargassi                                     
## Body_Length:genus_speciesChallengerosergia_hansjacobi                               
## Body_Length:genus_speciesChallengerosergia_talismani                                
## Body_Length:genus_speciesDeosergestes_corniculum                                    
## Body_Length:genus_speciesDeosergestes_henseni                                       
## Body_Length:genus_speciesGardinerosergia_splendens                                  
## Body_Length:genus_speciesParasergestes_armatus                                      
## Body_Length:genus_speciesParasergestes_vigilax                                      
## Body_Length:genus_speciesPhorcosergia_grandis                                       
## Body_Length:genus_speciesRobustasergia_robusta                                      
## Body_Length:genus_speciesRobustosergia_regalis                                      
## Body_Length:genus_speciesSergestes_atlanticus                                       
## Body_Length:genus_speciesSergia_tenuiremis                                          
## Body_Length:Preservationfresh                                                       
## Body_Length:Preservationparaformaldehyde                                            
## genus_speciesAllosergestes_sargassi:Preservationfresh                               
## genus_speciesChallengerosergia_hansjacobi:Preservationfresh                         
## genus_speciesChallengerosergia_talismani:Preservationfresh                          
## genus_speciesDeosergestes_corniculum:Preservationfresh                              
## genus_speciesDeosergestes_henseni:Preservationfresh                                 
## genus_speciesGardinerosergia_splendens:Preservationfresh                            
## genus_speciesParasergestes_armatus:Preservationfresh                                
## genus_speciesParasergestes_vigilax:Preservationfresh                                
## genus_speciesPhorcosergia_grandis:Preservationfresh                                 
## genus_speciesRobustasergia_robusta:Preservationfresh                                
## genus_speciesRobustosergia_regalis:Preservationfresh                                
## genus_speciesSergestes_atlanticus:Preservationfresh                                 
## genus_speciesSergia_tenuiremis:Preservationfresh                                    
## genus_speciesAllosergestes_sargassi:Preservationparaformaldehyde                    
## genus_speciesChallengerosergia_hansjacobi:Preservationparaformaldehyde              
## genus_speciesChallengerosergia_talismani:Preservationparaformaldehyde               
## genus_speciesDeosergestes_corniculum:Preservationparaformaldehyde                   
## genus_speciesDeosergestes_henseni:Preservationparaformaldehyde                      
## genus_speciesGardinerosergia_splendens:Preservationparaformaldehyde                 
## genus_speciesParasergestes_armatus:Preservationparaformaldehyde                     
## genus_speciesParasergestes_vigilax:Preservationparaformaldehyde                     
## genus_speciesPhorcosergia_grandis:Preservationparaformaldehyde                      
## genus_speciesRobustasergia_robusta:Preservationparaformaldehyde                     
## genus_speciesRobustosergia_regalis:Preservationparaformaldehyde                     
## genus_speciesSergestes_atlanticus:Preservationparaformaldehyde                      
## genus_speciesSergia_tenuiremis:Preservationparaformaldehyde                         
## Body_Length:genus_speciesAllosergestes_sargassi:Preservationfresh                   
## Body_Length:genus_speciesChallengerosergia_hansjacobi:Preservationfresh             
## Body_Length:genus_speciesChallengerosergia_talismani:Preservationfresh              
## Body_Length:genus_speciesDeosergestes_corniculum:Preservationfresh                  
## Body_Length:genus_speciesDeosergestes_henseni:Preservationfresh                     
## Body_Length:genus_speciesGardinerosergia_splendens:Preservationfresh                
## Body_Length:genus_speciesParasergestes_armatus:Preservationfresh                    
## Body_Length:genus_speciesParasergestes_vigilax:Preservationfresh                    
## Body_Length:genus_speciesPhorcosergia_grandis:Preservationfresh                     
## Body_Length:genus_speciesRobustasergia_robusta:Preservationfresh                    
## Body_Length:genus_speciesRobustosergia_regalis:Preservationfresh                    
## Body_Length:genus_speciesSergestes_atlanticus:Preservationfresh                     
## Body_Length:genus_speciesSergia_tenuiremis:Preservationfresh                        
## Body_Length:genus_speciesAllosergestes_sargassi:Preservationparaformaldehyde        
## Body_Length:genus_speciesChallengerosergia_hansjacobi:Preservationparaformaldehyde  
## Body_Length:genus_speciesChallengerosergia_talismani:Preservationparaformaldehyde   
## Body_Length:genus_speciesDeosergestes_corniculum:Preservationparaformaldehyde       
## Body_Length:genus_speciesDeosergestes_henseni:Preservationparaformaldehyde          
## Body_Length:genus_speciesGardinerosergia_splendens:Preservationparaformaldehyde     
## Body_Length:genus_speciesParasergestes_armatus:Preservationparaformaldehyde         
## Body_Length:genus_speciesParasergestes_vigilax:Preservationparaformaldehyde         
## Body_Length:genus_speciesPhorcosergia_grandis:Preservationparaformaldehyde          
## Body_Length:genus_speciesRobustasergia_robusta:Preservationparaformaldehyde         
## Body_Length:genus_speciesRobustosergia_regalis:Preservationparaformaldehyde         
## Body_Length:genus_speciesSergestes_atlanticus:Preservationparaformaldehyde          
## Body_Length:genus_speciesSergia_tenuiremis:Preservationparaformaldehyde             
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1457 on 381 degrees of freedom
## Multiple R-squared:  0.8953, Adjusted R-squared:  0.8774 
## F-statistic:  50.1 on 65 and 381 DF,  p-value: < 2.2e-16
#view anova table (see effects of variables)
anova(lm_pres)
## Analysis of Variance Table
## 
## Response: Eye_Diameter
##                                         Df Sum Sq Mean Sq   F value    Pr(>F)
## Body_Length                              1 55.549  55.549 2616.3984 < 2.2e-16
## genus_species                           13 10.906   0.839   39.5133 < 2.2e-16
## Preservation                             2  0.400   0.200    9.4204 0.0001016
## Body_Length:genus_species               13  1.256   0.097    4.5523 2.933e-07
## Body_Length:Preservation                 2  0.011   0.005    0.2496 0.7792470
## genus_species:Preservation              18  0.819   0.045    2.1418 0.0045283
## Body_Length:genus_species:Preservation  16  0.196   0.012    0.5765 0.9014187
## Residuals                              381  8.089   0.021                    
##                                           
## Body_Length                            ***
## genus_species                          ***
## Preservation                           ***
## Body_Length:genus_species              ***
## Body_Length:Preservation                  
## genus_species:Preservation             ** 
## Body_Length:genus_species:Preservation    
## Residuals                                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

The ANOVA table shows that, as we knew, body length, species, and their interaction all have significant effects on eye size (this is due to species specific patterns of ontogenetic scaling). However, preservation also has a significant effect (p=0.0001), as does the interaction between species and preservation (p = 0.005). This indicates that we should probably consider preservation error in our analyses, or that we may want to try to correct for any shrinkage in our data prior to analyses.

Phylogenetic distribution of allometric slopes

Below, we plot species-specific slopes for developmental allometry onto the sergestid tree to visualize the distribution and variation of the relationships between eye and body size throughout growth. There is potential to do more with these data; do we have predictions for what might drive species differences in developmental allometry? It looks like body size might be one factor, as the smaller-bodied species seem to have higher slopes. Whether species exhibit ontogenetic descent may be another factor to explore? We might predict that species that have ontogenetic descent and are going into dimmer habitats as they grow may invest in higher eye growth relative to body growth. Do all sergestids exhibit ontogenetic descent?

#Add slopes for developmental allometry into species dataframe ----

#add coefficient rownames to genus_species column
cc_species$genus_species <- rownames(cc_species)

#merge coefficients with species data by rowname
species_sma <- left_join(species, cc_species, by = "genus_species") %>%
  filter(genus_species != "Neosergestes_edwardsii")

# Prep data and phylogeny -----

#Make list of taxa to drop (in tree but not in dataset)
drops <- setdiff(tree$tip.label, as.character(species_sma$Binomial))

#Drop unwanted tips from phylogeny
tree.dev <- drop.tip(phy = tree, tip = drops) 

# Make Genus_species tip labels to replace Binomial codes
sp.tips <- data.frame(old = species_sma$Binomial, new = species_sma$genus_species)

# replace tip labels in phylogeny
tree.dev.plot <- sub.taxa.label(tree.dev, sp.tips)

#change row names of species data
rownames(species_sma) <- species_sma$genus_species

#check that phylogeny and data match exactly
name.check(tree.dev.plot, species_sma)

#resort trait dataset to the order of tree tip labels
species_sma <- species_sma[tree.dev.plot$tip.label, ] 

#make trait vector for slope
slope <- as.vector(species_sma$slope) 

#add tip label names to vector
names(slope) <- species_sma$genus_species

# Phylogeny with slopes for developmental eye-body allometry----
plotTree.wBars(tree.dev.plot, slope, 
               scale = 0.5, 
               tip.labels = TRUE,
               offset = 0.1,
               ftype = "bi",
               fsize = 1)

#export figure
png("../Figures/dev-phylo.png", width = 10, height = 7, units = "in", res = 500)

plotTree.wBars(tree.dev.plot, slope, 
               scale = 0.5, 
               tip.labels = TRUE,
               offset = 0.1,
               ftype = "bi",
               fsize = 1)

dev.off()

Evolutionary eye-body allometry in Sergestidae

We next examine interspecific eye-body allometry across sergestids using species means for eye diameter and body length. Here, we use phylogenetic least-squares regression (using the caper package) to fit the allometric relationship between eye size and body length across 15 species of sergestids.

Flagging to discuss calculation of species means with Lori

# Prep data ------

#check that names match in dataframe and tree
name.check(phy = tree, data = species, data.names = species$Binomial)
  
#set node names to null (or throws duplicate error)
tree$node.label <- NULL

#use caper function to combine phylogeny and data into one object (this function also matches species names in tree and dataset)
species.comp <- comparative.data(phy = tree, data = species, 
                            names.col = Binomial, vcv = TRUE, 
                            na.omit = FALSE, warn.dropped = TRUE)

#check for dropped tips or dropped species
species.comp$dropped$tips #phylogeny
species.comp$dropped$unmatched.rows #dataset

# PGLS model -------

#fit model
pgls_evo <- pgls(log10(Eye_diameter) ~ log10(Body_length), 
               data = species.comp, 
               lambda = "ML", #uses Maximum Liklihood estimate of lambda
               param.CI = 0.95)

#check model assumptions
par(mar = c(4,4,2,2))
par(mfrow = c(2, 2))
plot(pgls_evo)

par(mfrow = c(1, 1))

The PGLS regression shows that eyes scale with body length with a slope of 0.85, indicating that eyes have a negative/hypo allometry with body size. This is what is most common among all vertebrates and some invertebrates that have been examined so far (inverts tend to vary a lot more than verts).

The maximum likelihood estimate of lambda here is 1.00 (high phylogenetic signal), so that’s what caper has used to create the variance-covariance matrix for the regression. However, small sample sizes (< 20-30 data points) have no power to estimate lambda (see Freckleton et al. 2002), so tests for phylogenetic signal in this dataset are essentially meaningless (note that the 95% confidence interval for lambda here is 0 to 1, all possible values, and the p-values are non-signficant). For this reason, if any of our other comparative analyses estimate a lambda of 0, we should repeat them again forcing lambda to 1, just so we can see what difference that makes in the model fits (as we can’t know what lambda actually is in this small dataset). Alternatively, we could run all fits with lambda set to 0 (which ignores phylogeney) and again to 1 (highest phylogenetic signal) so that we can bin the possible results into the extremes.

#print model output 
summary(pgls_evo)
## 
## Call:
## pgls(formula = log10(Eye_diameter) ~ log10(Body_length), data = species.comp, 
##     lambda = "ML", param.CI = 0.95)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.03290 -0.01650  0.06354  0.07979  0.10039 
## 
## Branch length transformations:
## 
## kappa  [Fix]  : 1.000
## lambda [ ML]  : 1.000
##    lower bound : 0.000, p = 0.14833
##    upper bound : 1.000, p = 1    
##    95.0% CI   : (NA, NA)
## delta  [Fix]  : 1.000
## 
## Coefficients:
##                    Estimate Std. Error t value  Pr(>|t|)    
## (Intercept)        -1.21548    0.17275 -7.0358 8.855e-06 ***
## log10(Body_length)  0.77787    0.11327  6.8673 1.140e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.06792 on 13 degrees of freedom
## Multiple R-squared: 0.7839,  Adjusted R-squared: 0.7673 
## F-statistic: 47.16 on 1 and 13 DF,  p-value: 1.14e-05
#Likelihood plot for Pagel's lambda 
lambda.evo <- pgls.profile(pgls_evo, "lambda")
plot(lambda.evo, main = "Pagel's lambda", 
     sub = "Solid red line indicates estimate for lambda and broken red lines indcaite the 95% confidence interval")

#save coefficients of fits as object
cc_pgls_evo <- data.frame(coef(pgls_evo))

# Make plot ------
plot_evo <- ggplot(species, aes(x = Body_length, y = Eye_diameter, text = genus_species)) + 
  geom_point(size = 2, alpha = 0.6) + 
  scale_x_log10(name = "Body length (mm)") + #makes x axis log scale and named
  scale_y_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank()) +
  geom_abline(data = cc_pgls_evo, aes(intercept = cc_pgls_evo[1,1], slope = cc_pgls_evo[2,1])) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 1.5, nudge_x = 4, label = paste("slope =", round(cc_pgls_evo[2,1], digits = 2), sep = " "), size = 3) +
  geom_text(x = -Inf, y = Inf, hjust = -0.05, vjust = 3, label = paste("intercept =", round(cc_pgls_evo[1,1], digits = 2), sep = " "), size = 3)

#interactive plot
ggplotly(plot_evo)
# Export figure -----

Correlates of eye size in deep-sea sergestid shrimp

Next, we examine whether varition in relative eye size is correlated with ecological variables related to light regime and vision: depth range and bioluminescent organs.

Ecology and absolute and relative eye size across species

Here, we examine how absolute and relative eye size vary across the sergestid phylogeny, and whether these correlate to ecological variables. For visualizing relative eye size, we use the residuals of the PGLS fit for eye-body allometry, as these show eye size corrected by both body size and evolutionary relationships.

There are few enough species that it would be really nice to put an image of each species next to the phylogeny/associated data if you have them. If not photos, we could make vector images of species sillhouettes and scale them by body size for easy reference.

#Add residuals from PGLS fits to dataframe ----

#extract pgls residuals 
pglsres <- residuals(pgls_evo) 

#name residuals
colnames(pglsres) <- "pglsres" 

#merge residuals with species data by rowname
species.phy <- merge(species, pglsres, by = "row.names") 

#make column converted to observed eye size as factor of PGLS fit
species.phy$eyefactor <- 10^species.phy$pglsres

# Prep data and phylogeny -----

# Make Genus_species tip labels to replace Binomial codes
sp.tips <- data.frame(old = species.phy$Binomial, new = species.phy$genus_species)

# replace tip labels in phylogeny
tree.plot <- sub.taxa.label(tree, sp.tips)

#change row names of species data
rownames(species.phy) <- species.phy$genus_species

#check that phylogeny and data match exactly
name.check(tree.plot, species.phy)

#resort trait dataset to the order of tree tip labels
species.phy <- species.phy[tree.plot$tip.label, ] 

#make trait vector for absolute eye size
aveye <- as.vector(species.phy$Eye_diameter) 

#add tip label names to vector
names(aveye) <- species.phy$genus_species

#make trait vector for eye investment (PGLS residuals)
releye <- as.vector(species.phy$pglsres) #residuals of pgls
names(releye) <- species.phy$genus_species


# Phylogeny with eye size and investment -------

#color vector for photophores
col_ph <- c("lensed" = "palegreen",
            "none" = "honeydew",
            "pesta" = "cadetblue1",
            "unlensed" = "orchid1")

#plot tree with eye investment bars
plotTree.wBars(tree.plot, releye, 
               scale = 4, 
               tip.labels = TRUE,
               offset = 0.8,
               ftype = "bi",
               fsize = 1.2)

#add tip labels for absolute eye size
tiplabels(cex = 1.5*aveye,
          pch = 19, #shape of labels
          col = "mediumpurple",
          offset = 0) 

#add tip labels for eye size
tiplabels(text = species.phy$Organ,
          cex = 0.7,
          #frame = "none",
          bg = col_ph[species.phy$Organ],
          adj = 0,
          offset = 2) 
Phylogeny of the 15 sertestids sampled depicting absolute eye diameter (purple circles scaled by absolute size), relative eye investment (bars scaled by the residuals from the PGLS of eye diameter vs. body length), and photophore complexity (colored boxes).

Phylogeny of the 15 sertestids sampled depicting absolute eye diameter (purple circles scaled by absolute size), relative eye investment (bars scaled by the residuals from the PGLS of eye diameter vs. body length), and photophore complexity (colored boxes).

#Export figure------
png("../Figures/phylo.png", width = 8, height = 5.5, units = "in", res = 500)

#plot tree with eye investment bars
plotTree.wBars(tree.plot, releye, 
               scale = 4, 
               tip.labels = TRUE,
               offset = 0.8,
               ftype = "bi",
               fsize = 1.2)

#add tip labels for absolute eye size
tiplabels(cex = 1.5*aveye,
          pch = 19, #shape of labels
          col = "mediumpurple",
          offset = 0) 

#add tip labels for eye size
tiplabels(text = species.phy$Organ,
          cex = 0.7,
          #frame = "none",
          bg = col_ph[species.phy$Organ],
          adj = 0,
          offset = 2) 

dev.off()

Visualizing eye size and ecology

The following are exploratory plots that do not account for the phylogenetic distribution of traits. Plots are interactive: hover over a data point to see the species and x/y values.

Photophore complexity

Prediction: We might predict that species that have highly developed photophrores may be more likely to engage in visual signaling, and thus may have higher investments in eye size. However, this is the opposite of what we see here (species with organs of pesta have the smallest eyes and the lowest relative investment in eye size). As BL can also be used for camouflage, and BL from other species is ubiquitous, there is a reasonable expectation that they could be uncorrelated as well.

#Boxplots for absolute eye size across photophores -----
plot_abs_photo <- ggplot(data = species.phy, aes(x = reorder(Organ, Eye_diameter, FUN = mean), y = Eye_diameter)) + 
  geom_boxplot(notch = FALSE, outlier.shape = NA, alpha = 0.6) + #controls boxes
  stat_summary(fun.y = mean, colour = "black", geom = "point", shape = 18, size = 3, show_guide = FALSE) + #controls what stats shown
 geom_jitter(shape = 19, size = 1.3, alpha = 0.5, position = position_jitter(0.15)) +
  theme(text = element_text(size=14), panel.background = element_blank(), axis.line = element_line(colour = "black")) + #controls background
  xlab("Photophore complexity") +
  ylab("Eye diameter (mm)") +
  ggtitle("Absolute eye size")

ggplotly(plot_abs_photo)
#Boxplots for relative eye investment across photophores -----
plot_rel_photo <- ggplot(data = species.phy, aes(x = reorder(Organ, pglsres, FUN = mean), y = pglsres)) + 
  geom_boxplot(notch = FALSE, outlier.shape = NA, alpha = 0.6) + #controls boxes
  stat_summary(fun.y = mean, colour = "black", geom = "point", shape = 18, size = 3, show_guide = FALSE) + #controls what stats shown
 geom_jitter(shape = 19, size = 1.3, alpha = 0.5, position = position_jitter(0.15)) +
  theme(text = element_text(size=14), panel.background = element_blank(), axis.line = element_line(colour = "black")) + #controls background
  xlab("Photophore complexity") +
  ylab("Relative eye investment") +
  ggtitle("Relative eye investment")

ggplotly(plot_rel_photo) 

Maximum depth

Prediction: We might expect that eye size would increase with maximum depth to some point, then may level off and become uncorrelated once the pelagic zone occupied lacks downwelling light for vision. Maximum recorded depth isn’t the best indicator of where an animal actually spends its time, so we could expect that the actual depth occupied at the brightest time of day could be somewhere a bit shallower than these values, so tough to say whether and where leveling off might occur.

Visual inspection of the data shows that absolute eye size does seem to increase with maximum depth. However, body size also increases with maximum depth, and both relative eye size (eye/body) and eye investment (residuals of PGLS for eye ~ body) show litte to no correlation with maximum depth.

# Absolute eye size by max depth ------
plot_abs_max <- ggplot(species.phy, aes(x = Eye_diameter, y = Max, text = genus_species)) + 
  geom_point(size = 2, alpha = 0.6) + 
  scale_y_reverse(name = "Maximum depth (m)") + #makes x axis log scale and named
  scale_x_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())

ggplotly(plot_abs_max)
# Body size with max depth -----
plot_body_max <- ggplot(species.phy, aes(x = Body_length, y = Max, text = genus_species)) + 
  geom_point(size = 2, alpha = 0.6) + 
  scale_y_reverse(name = "Maximum depth (m)") + #makes x axis log scale and named
  scale_x_log10(name = "Body length (mm)") + #makes y axis log scale and named
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())

ggplotly(plot_body_max)
# Eye size/body size with max depth -----
plot_eb_max <- ggplot(species.phy, aes(x = Relative_size, y = Max, text = genus_species)) + 
  geom_point(size = 2, alpha = 0.6) + 
  scale_y_reverse(name = "Maximum depth (m)") + #makes x axis log scale and named
  scale_x_log10(name = "Eye:body") + #makes y axis log scale and named
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())

ggplotly(plot_eb_max)
# Relative eye investment by max depth -----
plot_inv_max <- ggplot(species.phy, aes(x = eyefactor, y = Max, text = genus_species)) + 
  geom_point(size = 2, alpha = 0.6) + 
  scale_y_reverse(name = "Maximum depth (m)") + #makes x axis log scale and named
  scale_x_log10(name = "Relative eye investment") + #makes y axis log scale and named
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())

ggplotly(plot_inv_max)
# Absolute eye size by depth range -----
plot_abs_range <- ggplot(species.phy, aes(x = Eye_diameter, y = (Max - Min), text = genus_species)) + 
  geom_point(size = 2, alpha = 0.6) + 
  scale_y_reverse(name = "Depth range (m)") + #makes x axis log scale and named
  scale_x_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())

ggplotly(plot_abs_range)
# Relative eye investment by depth range -----
plot_inv_range <- ggplot(species.phy, aes(x = eyefactor, y = (Max - Min), text = genus_species)) + 
  geom_point(size = 2, alpha = 0.6) + 
  scale_y_reverse(name = "Depth range (m)") + #makes x axis log scale and named
  scale_x_log10(name = "Relative eye investment") + #makes y axis log scale and named
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())

ggplotly(plot_inv_range)
# Absolute eye size by mid depth -----

#add column for midpoint of depth range
species.phy <- species.phy %>%
  mutate(depth_mid = ((Max - Min)/2) + Min)

#plot
plot_abs_mid <- ggplot(species.phy, aes(x = Eye_diameter, y = depth_mid, text = genus_species)) + 
  geom_point(size = 2, alpha = 0.6) + 
  scale_y_reverse(name = "Midpoint depth (m)") + #makes x axis log scale and named
  scale_x_log10(name = "Eye diameter (mm)") + #makes y axis log scale and named
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())

ggplotly(plot_abs_mid)
# Relative eye investment by mid depth -----
plot_inv_mid <- ggplot(species.phy, aes(x = eyefactor, y = depth_mid, text = genus_species)) + 
  geom_point(size = 2, alpha = 0.6) + 
  scale_y_reverse(name = "Midpoint depth (m)") + #makes x axis log scale and named
  scale_x_log10(name = "Relative eye investment") + #makes y axis log scale and named
  theme_bw() +
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())

ggplotly(plot_inv_mid)

Comparative models

Eye diameter ~ body length * photophore complexity

Below, I run a PGLS model of log10(eye diameter) vs. log10(body length) with photophore complexity as a covariate (= phylogenetic ANCOVA). Because we have so many states for photophore complexity and a small sample of species, we cannot look at the interaction term for body size x photophore complexity because then our states outnumber our number of observations.

#run PGLS modelusing the Maximum Liklihood estimate of lambda
pgls_photo <- pgls(log10(Eye_diameter) ~ log10(Body_length) + Organ, 
               data = species.comp, 
               lambda = "ML",
               #uses Maximum Liklihood estimate of lambda
               param.CI = 0.95)

#evaluate model assumptions
par(mfrow = c(2,2)) #makes your plot window into 2x2 panels
plot(pgls_photo) #plot the linear model

par(mfrow = c(1,1)) #set plotting window back to one plot

Model outputs

This output from anova(model) gives the overall significance of the main, 2-way, 3-way, etc. interactions and sequential sums of squares. Here we can see that body length has a significant effect on eye size (expected, and we already knew that), but we also find that photophore complexity has a significant effect on eye size when the effect of body size is accounted for (p = 0.01).

#anova
anova(pgls_photo)
## Analysis of Variance Table
## Sequential SS for pgls: lambda = 0.00, delta = 1.00, kappa = 1.00
## 
## Response: log10(Eye_diameter)
##                    Df  Sum Sq Mean Sq  F value    Pr(>F)    
## log10(Body_length)  1 0.33032 0.33032 160.2202 1.766e-07 ***
## Organ               3 0.03640 0.01213   5.8858   0.01397 *  
## Residuals          10 0.02062 0.00206                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

This output from summary(model) gives PGLS coefficients as well significance with respect to contrasts. Note that for a categorical variable, R will treat whichever factor comes first as the “control” and will test other factors for significant differences relative to that “control”. Here, the comparison state (“control”) for photophore complexity is set to “unlensed” (estimate listed under (Intercept) - that’s why you don’t see that state as a comparison line) and the other states show the difference in the estimated intercept for that state compared to the unlensed state, and whether that difference is significant. You can see that “none” and “lensed” show no significant difference in intercept from “unlensed”, but that “pesta” does show a significant difference (p = 0.01) and has a lower intercept, indicating smaller relative eye size.

The output also includes the maximum-likelihood estimation of lambda, which represents the phylogenetic signal for PGLS model residuals. Here, the estimate for lambda is 0.00, which means that there is no phylogenetic signal. This doesn’t mean we shouldn’t include phylogeny in the model (no reason not to), but it does mean that it’s the same result we would get if we assessed this without including phylogeny. This also doesn’t mean that there is not phylogenetic signal in the individual traits that we put into the model; this is only lambda for the model residuals, not the traits themselves.

#summary
summary(pgls_photo)
## 
## Call:
## pgls(formula = log10(Eye_diameter) ~ log10(Body_length) + Organ, 
##     data = species.comp, lambda = "ML", param.CI = 0.95)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.068782 -0.021392 -0.003054  0.028138  0.038374 
## 
## Branch length transformations:
## 
## kappa  [Fix]  : 1.000
## lambda [ ML]  : 0.000
##    lower bound : 0.000, p = 1    
##    upper bound : 1.000, p = 0.10592
##    95.0% CI   : (NA, NA)
## delta  [Fix]  : 1.000
## 
## Coefficients:
##                     Estimate Std. Error t value  Pr(>|t|)    
## (Intercept)        -0.940877   0.152641 -6.1640 0.0001064 ***
## log10(Body_length)  0.647943   0.096499  6.7145 5.269e-05 ***
## Organnone          -0.032417   0.056619 -0.5726 0.5795934    
## Organpesta         -0.115872   0.037867 -3.0600 0.0120440 *  
## Organunlensed       0.025402   0.041107  0.6180 0.5504201    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.04541 on 10 degrees of freedom
## Multiple R-squared: 0.9468,  Adjusted R-squared: 0.9255 
## F-statistic: 44.47 on 4 and 10 DF,  p-value: 2.449e-06

Eye diameter ~ body length + maximum depth

Below, I run a PGLS model of klg10(eye diameter) vs. log10(body length) with maximum depth as a covariate. This model does not test for the interaction between body length and maximum depth.

Question to consider: Should we log maximum depth in these analyses? Allometric data are typically logged and I have done so throughout these analyses. Depth you could argue either way, but since light attenuates exponentially with depth and we are thinking about light environment here I’d be inclined to log? Currently, eye size and body size are logged and depth is unlogged in the below analysis.

#run PGLS modelusing the Maximum Liklihood estimate of lambda
pgls_depth <- pgls(log10(Eye_diameter) ~ log10(Body_length) + Max, 
               data = species.comp, 
               lambda = "ML", bounds = list(lambda = c(1e-05, 1)),
               #uses Maximum Liklihood estimate of lambda
               param.CI = 0.95)

#evaluate model assumptions
par(mfrow = c(2,2)) #makes your plot window into 2x2 panels
plot(pgls_depth) #plot the linear model

par(mfrow = c(1,1)) #set plotting window back to one plot

Model outputs

This output from anova(model) gives the overall significance of the main effects and sequential sums of squares. Here we can see that body length has a significant effect on eye size when depth is corrected for (expected, and we already knew that), but that maximum depth also has a significant effect on eye size corrected for the effects of body size (p < 0.05).

Note: I originally ran this model with an additional interaction term for depth x body length, and there was not a significant effect for maximum depth (p = 0.05) or the interaction between max depth and body length (p = 0.67). Note though that depth was very close to being significant, so this may result from our low species sample size (n = 15) combined with additional model paramaterization yielding less statistical power; additional species sampling may yield a different result.

Second note: There is a postitive relationship between body size and maximum depth in these data (body length increases with maximum depth). This raises some concern about multicolinearity, as these models assume that our two predictive variables are independent of each other (and in this case, that is at least partially false). The finding that species mean body size increases with maximum depth could be a potentially interesting thing to have a look at as well, there are a number of studies looking at body size with depth in the deep ocean so it’s definitely of interest, I think especially people thinking about carbon/food availability.

#anova
anova(pgls_depth)
## Analysis of Variance Table
## Sequential SS for pgls: lambda = 0.00, delta = 1.00, kappa = 1.00
## 
## Response: log10(Eye_diameter)
##                    Df  Sum Sq Mean Sq F value    Pr(>F)    
## log10(Body_length)  1 0.33032 0.33032 147.601 4.214e-08 ***
## Max                 1 0.03017 0.03017  13.479  0.003199 ** 
## Residuals          12 0.02685 0.00224                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

This output from summary(model) gives PGLS coefficients as well significance with respect to contrasts.

The output also gives a maximum-likelihood estimation of lambda = 0.00, indicating that there is no phylogenetic signal in the residuals of the model (and that caper used a lambda of 0 to set the variance-covariance matrix in the regression). However, small sample sizes have no power to estimate the true value of lambda, so this is not very meaningful in this dataset. We should decide on a strategy for analysis (as mentioned above, may want to run all models again with lambda fixed at 1 as a comparison for the supp info).

#summary
summary(pgls_depth)
## 
## Call:
## pgls(formula = log10(Eye_diameter) ~ log10(Body_length) + Max, 
##     data = species.comp, lambda = "ML", param.CI = 0.95, bounds = list(lambda = c(1e-05, 
##         1)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.074077 -0.025220  0.001858  0.030619  0.070874 
## 
## Branch length transformations:
## 
## kappa  [Fix]  : 1.000
## lambda [ ML]  : 0.000
##    lower bound : 0.000, p = 1    
##    upper bound : 1.000, p = 0.3552
##    95.0% CI   : (NA, NA)
## delta  [Fix]  : 1.000
## 
## Coefficients:
##                       Estimate  Std. Error t value  Pr(>|t|)    
## (Intercept)        -1.1603e+00  1.2761e-01 -9.0929 9.905e-07 ***
## log10(Body_length)  6.6745e-01  9.6945e-02  6.8849 1.687e-05 ***
## Max                 9.8318e-05  2.6779e-05  3.6714  0.003199 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.04731 on 12 degrees of freedom
## Multiple R-squared: 0.9307,  Adjusted R-squared: 0.9191 
## F-statistic: 80.54 on 2 and 12 DF,  p-value: 1.111e-07
#Here, the fit for predicting log10(eye size) is -1.18 + log10(0.64 * Body length) + log10(0.0001 * Max depth)

Eye diameter ~ body length * depth range

Below, I run a PGLS model of log10(eye diameter) vs. log10(body length) with depth range as a covariate. This model does not test for the interaction between body length and depth range.

#run PGLS modelusing the Maximum Liklihood estimate of lambda
pgls_range <- pgls(log10(Eye_diameter) ~ log10(Body_length) + (Max - Min), 
               data = species.comp, 
               lambda = "ML", bounds = list(lambda = c(1e-05, 1)),
               #uses Maximum Liklihood estimate of lambda
               param.CI = 0.95)

#evaluate model assumptions
par(mfrow = c(2,2)) #makes your plot window into 2x2 panels
plot(pgls_range) #plot the linear model

par(mfrow = c(1,1)) #set plotting window back to one plot

Model outputs

This output from anova(model) gives the overall significance of the main effects and sequential sums of squares. Both body length (as expected) and depth range (p < 0.05) show significant effects on eye size.

#anova
anova(pgls_range)
## Analysis of Variance Table
## Sequential SS for pgls: lambda = 0.00, delta = 1.00, kappa = 1.00
## 
## Response: log10(Eye_diameter)
##                    Df  Sum Sq Mean Sq F value    Pr(>F)    
## log10(Body_length)  1 0.33032 0.33032 147.601 4.214e-08 ***
## Max                 1 0.03017 0.03017  13.479  0.003199 ** 
## Residuals          12 0.02685 0.00224                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

This output from summary(model) gives PGLS coefficients as well significance with respect to contrasts.

The output also gives a maximum-likelihood estimation of lambda = 0.00, indicating that there is no phylogenetic signal in the residuals of the model (and that caper used a lambda of 0 to set the variance-covariance matrix in the regression). However, small sample sizes have no power to estimate the true value of lambda, so this is not very meaningful in this dataset. We should decide on a strategy for analysis (as mentioned above, may want to run all models again with lambda fixed at 1 as a comparison for the supp info).

You can see below that the coefficients are really similar to those for max depth. I think that’s because in many cases depth range = max depth (min depth only varies between 0 and 300 m). They are both probably showing pretty much the same thing here.

#summary
summary(pgls_range)
## 
## Call:
## pgls(formula = log10(Eye_diameter) ~ log10(Body_length) + (Max - 
##     Min), data = species.comp, lambda = "ML", param.CI = 0.95, 
##     bounds = list(lambda = c(1e-05, 1)))
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.074077 -0.025220  0.001858  0.030619  0.070874 
## 
## Branch length transformations:
## 
## kappa  [Fix]  : 1.000
## lambda [ ML]  : 0.000
##    lower bound : 0.000, p = 1    
##    upper bound : 1.000, p = 0.3552
##    95.0% CI   : (NA, NA)
## delta  [Fix]  : 1.000
## 
## Coefficients:
##                       Estimate  Std. Error t value  Pr(>|t|)    
## (Intercept)        -1.1603e+00  1.2761e-01 -9.0929 9.905e-07 ***
## log10(Body_length)  6.6745e-01  9.6945e-02  6.8849 1.687e-05 ***
## Max                 9.8318e-05  2.6779e-05  3.6714  0.003199 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.04731 on 12 degrees of freedom
## Multiple R-squared: 0.9307,  Adjusted R-squared: 0.9191 
## F-statistic: 80.54 on 2 and 12 DF,  p-value: 1.111e-07
#Here, the fit for predicting log10(eye size) is -1.18 + log10(0.64 * Body length) + log10(0.0001 * Max depth)

Thoughts on analyses so far

It looks like we are detecting effects of both depth and photophore complexity on relative eye size, which is super cool!! In terms of other models to build, since we have a low number of species for comparative methods (n = 15), I’m inclined to avoid building larger additive models for these effects, as that will greatly increase the ratio of parameter number to sample size. That’s just my thought though, we could certainly try it if you want to, as well as look at alternate additive approaches like ranking models by AIC. With PGLS, you have to fix lambda across the models instead of estimating via max likelihood or else the models aren’t comparable via AIC, but since the estimate of lambda is unreliable anyway here we could try repeating analyses with lambda fixed at 0 and 1 to bin the range of possible results depending how much phylogenetic signal is truly present.

I think the best approach is to use these models to do specific hypothesis testing. We can certainly plug lots of variables into the models in various combinations and do model comparisons to find the ‘best fits’, but the more things we test, the more likely we are to find things that come up significant, and the more parameters we add, the more likely higher-parameter models will fit better. So, in general I’d advocate for only testing variables that you have a hypothesis for in terms of why that might affect relative eye size, or things that you see interesting patterns in when looking at the data.

I’m super excited about your data and initial results and am always happy to add plots/analyses here in R as you think of more cool things to look at!!